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Abstract 

This  program  concerned  the  development  and  analysis  of  nanoenergetic  and  high  hydrogen 
content  materials  for  use  in  propellants.  During  the  first  two  years  of  this  program,  research  was 
conducted  on:  1)  experimental  and  multi-scale  modeling  of  the  thermo-chemical  behavior  of  nano¬ 
aluminum  particles,  2)  combustion  and  decomposition  studies  of  metal  hydrides  including  ammonia 
borane,  and  3)  the  integration  of  nanomaterials  into  propellants  and  the  combustion  analysis  of  these 
materials.  During  the  third  and  final  year  of  the  program,  the  research  concerned  the  development  of 
theoretical  methods  for  understanding  and  predicting  anisotropic  thermal  transport  and  energy  release  in 
rocket  propellant  formulations.  This  report  describes  accomplishments  from  the  effort. 
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I.  Introduction 

The  last  few  decades  have  witnessed  a  tremendous  increase  in  space-related  activities  and  many 
more  challenging  space  missions  are  currently  being  planned.  So  far,  most  space  exploration  and 
utilization  activities  have  primarily  employed  chemical  propulsion  technology  and  the  situation  will 
continue  in  the  foreseeable  future.  The  progressively  ambitious  future  space-mission  goals  have  identified 
an  acute  need  for  novel,  better  performing,  cost-effective,  and  safer  ways  of  space  exploration  and 
utilization.  Unfortunately,  none  of  the  current  operational  propellants  have  all  the  desired  properties,  and 
the  selection  of  a  particular  propellant  combination  is  largely  a  compromise  based  on  various  factors  such 
as  physical  properties,  ignition  and  combustion  characteristics,  cost,  safety,  environmental  hazards, 
processing,  material  compatibility,  and  availability  [1].  Substantial  efforts  are  thus  needed  towards  the 
development  of  new  propellants  to  meet  the  requirements  of  the  next  generation  of  space  exploration  and 
utilization  activities.  The  focus  on  advancing  the  existing  propellant  technologies,  at  this  juncture,  will 
benefit  all  major  space  programs  for  the  Air  Force,  NASA,  and  industry.  A  typical  space  mission  requires 
two  different  types  of  propulsion  systems,  viz.,  launch  vehicle  (i.e.  primary  and  upper  stages)  and 
spacecraft  (i.e.  on-board)  systems.  The  combination  of  solid  and  liquid  rockets  is  typical  of  launch 
systems,  wherein  solid  propellants  provide  a  cost  effective  and  simple  way  to  deliver  very  high  thrust 
required  at  lift-off,  and  liquid  propellants  fulfill  the  higher  specific-impulse  requirement  of  upper  stages. 
Currently,  small  launches  use  solid  propellants  for  the  sake  of  economy,  simplicity,  and  operability.  The 
medium  and  heavy  launches  use  liquid  propellants  in  the  main  propulsion  system,  and  are  usually  assisted 
by  solid  rocket  boosters  in  the  initial  stage  to  meet  the  thrust  requirements. 

The  propulsion  system  (including  propellants),  which  is  often  the  largest  and  the  most  massive 
component  of  a  launch  vehicle,  provides  a  powerful  leverage  that  can  be  gained  by  using  high  density  and 
better-performing  propellants.  The  overall  objective  for  the  advancement  of  such  systems  is  crucially 
dependent  on  the  density  and  specific  impulse  of  the  propellant.  Conventionally,  the  requirement  of 
higher  specific  impulse  has  led  to  the  use  of  cryogenic  bipropellants. 

Liquid  oxygen  (LOX)  and  liquid  hydrogen  (LH)  have  been  most  widely  used  for  space  missions 
[2],  especially  in  the  western  world.  Some  of  the  advantages  in  using  LOX/LH  propellants  are: 

•  low  molecular  weight  of  hydrogen  which  implies  higher  specific  impulse  (30%  -  40%  higher 

than  other  existing  rocket  propellants) 

•  non-toxic  and  clean  combustion  exhaust 

•  superior  chamber  cooling  capability 

•  relatively  stable  combustion 

On  the  other  hand,  the  disadvantages  associated  with  the  use  of  LOX/LH  cryogenic  bipropellants 
are: 

•  low  density  of  hydrogen  which  necessitates  a  large  storage  volume 

•  bulky  tanks  leading  to  increase  in  the  structural  mass 

•  safety  issues  linked  with  hydrogen  and  oxygen 

•  very  low  storage  temperatures 

•  propellant  compatibility  with  materials  of  storage  containers 

Due  to  low  temperatures  of  cryogenic  propellants,  it  is  difficult  to  store  them  for  long  periods  of 
time.  Consequently,  such  propellants  are  less  desirable  for  use  in  rockets  that  must  be  kept  launch  ready 
for  months  together.  The  tank  for  storing  liquid  propellants  is  typically  the  highest  nonexpendable  mass  in 
a  chemical  propulsion  system.  Since  liquid  hydrogen  has  a  very  low  density,  the  propellant  volume 
required  is  many  times  greater  than  other  denser  fuels,  thus  requiring  bulky  and  massive  storage  tanks. 
The  use  of  high-density  propellants  is  an  excellent  way  to  reduce  the  structural  mass  of  storage  tanks  and 
associated  propellant-conditioning  and  flow-control  components.  The  increase  in  propellant  density  thus 
translates  into  a  lower  take-off  weight  and  a  larger  vehicle  payload  capacity.  Higher-energy  propellant 
delivers  better  specific  impulse,  and  consequently  helps  reduce  the  required  propellant  weight  for  a  given 
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total  impulse.  Such  propellants  can  achieve  the  desired  performance  of  a  mission  with  smaller  and  lighter 
launch  vehicles. 

There  exists  a  need  to  thoroughly  identify,  and  investigate  novel  propellants  in  order  to  achieve 
desired  combinations  with  higher  density,  higher  specific  impulse,  and  better  thermal  properties. 
Alternatives  to  such  propellant  systems  as  LOX/LH  need  to  be  explored  in  a  systematic  manner. 

Additives  to  hydrogen  or  the  use  of  hydrocarbons  may  allow  upper  stages  to  deliver  a  better 
overall  performance  with  reduced  volume  and  cost  [2].  One  such  fuel  with  desirable  properties  which 
may  replace  LH  is  liquid  methane.  In  the  past  forty  years,  LOX/methane  is  the  only  new  propellant 
combination  that  has  been  adopted  for  flight  engines  in  the  United  States.  Such  a  lack  of  innovation 
highlights  the  fact  that  the  advancements  in  propellant  technology  have  been  found  wanting  for  a  long 
time.  Although  LOX/methane  provides  a  reasonable  propulsive  performance  for  certain  applications, 
more  appropriate  and  effective  propellant  combinations  must  be  continuously  explored,  especially  in  light 
of  the  encouraging  advances  in  the  field  of  energetic  materials  in  the  past  twenty  years.  Ammonia  borane 
has  already  been  used  with  success  to  produce  hydrogen  for  chemical  lasers  and  fuel  cells.  The 
alternatives  being  considered  for  the  oxidizer  LOX  are  hydroxylammonium  nitrate  (HAN),  ammonium 
dinitramide  (ADN),  hydrazinium  nitroformate  (HNF),  and  water  (H20).  Gelled  propellants  are  being 
tested  to  determine  the  practicality  of  letting  propellants  freeze  at  low  environmental  temperatures  and 
thawing  them  only  when  required  for  use.  In  addition  to  the  above  possible  combinations,  there  is  a 
concept  to  harness  powerful  chemical  bonds  between  individual  atoms  of  hydrogen,  boron,  carbon,  and 
aluminum.  The  atoms  could  be  arrested  with  a  cryogenic  solid,  and  released  as  they  enter  the  rocket 
engine.  The  ongoing  progress  in  nanotechnology  can  make  such  propellants  a  reality. 

Among  the  commonly  used  bi-propellants  are  monomethylhydrazine  (CH3NHNH2, 
MMH)/nitrogen  tetroxide  (N204,  NTO)  and  hydrazine  (N2H4)/NTO.  Due  to  the  toxicity  of  these 
propellants,  other  bi-propellants  such  as  LOX/ethanol  and  H202/hydrocarbon  are  being  explored.  The 
catalytic-decomposed  N2H4  has  been  used  extensively  as  a  monopropellant  for  the  primary  propulsion  of 
smaller  spacecrafts  as  well  as  for  gas-generator  applications  [3].  H202  has  also  been  used  as  a 
monopropellant  in  the  past.  With  a  goal  to  develop  monopropellants  with  significantly  better  performance 
and  no  hazardous  properties,  efforts  have  been  focused  on  aqueous  solution  of  HAN.  The  HAN-based 
monopropellants  have  higher  densities  and  lower  freezing  points  than  the  state-of-the-art  hydrazine. 
Furthermore  they  pose  no  vapor  hazard,  and  do  not  require  any  extraordinary  storage,  handling,  or 
disposal  procedures.  For  very  small  and  micro  spacecrafts,  several  alternative  propulsion  technologies  can 
provide  performance  and  system  benefits.  Examples  include  (a)  a  warm  gas  (mixture  of  hydrogen, 
oxygen,  and  an  inert  gas)  propulsion  system  that  offers  a  high  specific  impulse  as  compared  to  the 
conventional  cold-gas  systems,  (b)  exothermic  decomposition  of  solid  and  hybrid  systems,  which  offer 
the  high  density  and  simplicity  of  solid  propellants  for  low-thrust,  and  quick-response  applications,  and 
(c)  bipropellant  systems  using  microelectronic  mechanical  system  (MEMS)  fabrication  technology.  In  the 
past,  the  stringent  system  design  requirements  and  the  reluctance  to  change  them  have  obviated  the  use 
and  application  of  promising  high-energetic  materials.  But  it  has  become  increasingly  clear  that  to  meet 
the  needs  of  future  space-propulsion  activities,  a  stage  has  been  reached  wherein  novel  and  well 
characterized  HDEMs  may  dictate  the  overall  propulsion  system  design. 

This  research  effort  investigated  new  energetic  materials  that  feature  nanometer  length  scales 
and/or  materials  with  potential  of  yielding  high  hydrogen  yields.  In  a  previous  effort,  frozen  propellants 
based  on  nAl  and  ice  were  studied,  which  because  of  the  nanometer  sized  aluminum  particles  resulted  in 
fast  burning  rates  and  produced  hydrogen  gas  as  the  primary  gaseous  combustion  product.  The  present 
program  continued  research  efforts  on  nanometer  metallics  as  additives  to  propellants,  with  particular 
emphasis  on  nanometer  aluminum,  and  also  considered  materials  with  high  hydrogen  content,  such  as 
ammonia  borane  and  several  metal  hydrides.  During  the  first  two  years  of  this  effort,  research  was 
conducted  at  Penn  State  University,  Purdue  University,  and  Georgia  Tech.  These  efforts  included:  1) 
multi-scale  modeling  of  the  thermo-chemical  behavior  of  nano-aluminum  particles,  2)  further 
development  of  aluminum  and  water  solid  propellants,  3)  performance  analysis  of  metal  hydride  fuel 
additives,  4)  encapsulation  of  nanoscale  particles  in  oxidizer  crystals,  5)  investigation  of  inclusion 
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materials  in  micron  scale  aluminum,  6)  combustion  analysis  of  nano-aluminum  with  hydroxyl  ammonium 
nitrate  (HAN)/water  mixtures,  7)  decomposition  and  oxidation  studies  of  ammonia  borane  (with 
Professors  Adri  van  Duin  and  Stefan  Thynell  at  Penn  State),  and  8)  studies  on  the  usage  of  several  of 
these  materials  in  hybrid  rocket  motors.  Examples  of  this  research  are  presented  in  the  sections  that 
follow  with  references  to  other  papers  published  in  the  literature  as  a  result  of  this  effort.  During  the  third 
and  final  year  of  this  program,  the  research  was  conducted  at  the  University  of  Missouri-Columbia, 
Caltech,  and  the  University  of  Illinois  Urbana-Champaign.  This  research  focused  on  the  development  of 
theoretical  methods  for  understanding  and  predicting  anisotropic  thermal  transport  and  energy  release  in 
rocket  propellant  formulations.  A  summary  of  the  accomplishments  of  this  one  year  effort  are  also 
provided  in  the  report. 
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II.  Ammomia  Borane  as  a  High  Hydrogen  Content  Propellant  Additive  (Penn  State) 

Future  space  missions  with  chemical  propulsion  systems  will  require  the  development  of 
propellants  with  increased  performance.  To  increase  the  energy  density  (per  unit  mass)  over  current 
hydrocarbons,  earlier  studies  identified  hydrogen,  beryllium,  boron,  or  some  combination  thereof  as 
elements  of  future  fuels  [4].  Molecular  hydrogen  has  been  used  successfully  as  a  liquid  fuel  for  rocket 
propulsion.  However,  it  requires  cryogenic  facilities,  which  complicates  handling  and  rocket  design. 
Beryllium  and  its  compounds  have  generally  been  deemed  too  toxic  to  be  viable  fuels.  Boron,  however, 
is  of  strong  interest  due  its  high  combustion  enthalpy  on  both  a  mass  and  volume  basis.  It  is  typically 
introduced  as  a  particle  with  a  passivation  layer,  and  as  such  there  are  significant  difficulties  with  ignition 
and  complete  combustion  [5].  Similarly,  boranes  (BxHy)  have  been  of  interest  for  years  as  a  propellant 
fuel  because  of  their  high  energy  density.  The  liquid  pentaborane  (B5H9)  has  60%  greater  energy  density 
than  conventional  jet  fuel  [4].  However,  problems  with  the  pyrophoricity  and  toxicity  of  boranes 
ultimately  proved  them  too  dangerous  for  widespread  use. 

Recently,  ammonia  borane  (AB,  NH3BH3),  a  solid  compound  at  room  temperature,  has  been  the 
subject  of  vigorous  study  to  examine  its  potential  for  storing  hydrogen  at  high  densities.  AB  contains 
19.6%  hydrogen  by  mass.  The  initial  studies  on  AB  decomposition  [6,  7]  showed  that  it  would  slowly 
decompose  to  produce  hydrogen  gas.  Unlike  the  previously  investigated  borane  fuels  (B2H6,  B5H9),  AB  is 
not  known  to  have  significant  toxicity  or  safety  problems,  and  it  does  not  require  any  special  storage 
conditions  [8]. 

Increasing  specific  impulse  ( Isp )  is  one  of  the  main  objectives  in  developing  novel  rocket 
propellants  for  deep  space  travel.  Since  Isp  is  proportional  to  the  square  root  of  chamber  temperature 
divided  by  the  average  molecular  weight  of  the  products,  AB  is  of  interest  because  it  undergoes  an 
exothermic  decomposition  to  release  H2  [9].  Thermochemical  calculations  and  preliminary  laboratory 
scale  motor  tests  confirm  that  AB  can  indeed  improve  specific  impulse  in  a  paraffin-fueled  hybrid  rocket 
motor  [10].  However,  this  study  also  showed  that  significant  additions  of  AB  result  in  a  considerable 
quantity  of  condensed-phase  products  accumulating  on  the  solid-fuel  surface,  leading  to  a  decreased 
regression  rate.  To  understand  the  chemistry  related  to  these  results,  a  better  understanding  of  AB 
decomposition  is  required. 

A.  Using  Molecular  Dynamics  Simulations  with  a  ReaxFF  Reactive  Force  Field  to 

Develop  a  Kinetic  Mechanism  for  Ammonia  Borane  Oxidation  (with  Adri  van  Duin,  Penn 

State  University) 

The  use  of  chemical  kinetic  modeling  has  been  used  with  great  success  in  understanding 
combustion  phenomena,  but  such  modeling  requires  the  knowledge  of  all  relevant  elementary  chemical 
reactions  and  the  reaction  rates  associated  with  them.  The  kinetics  of  boron  reactions  in  the  presence  of 
hydrocarbon  oxidation  [11-14]  have  been  reported  previously.  Conversely,  elementary  reactions  and  rate 
parameters  describing  the  oxidation  of  B/N/H  species  are  significantly  lacking  from  the  literature. 
Acquiring  experimental  data  for  elementary  reaction  rates  in  a  completely  novel  system  is  time 
consuming  and  often,  although  still  very  important.  Computational  techniques,  such  as  molecular 
dynamics  (MD),  have  matured  along  with  modern  computing  power  and  can  aid  greatly  in  the 
development  of  kinetic  mechanisms  for  novel  systems. 

Advances  in  the  speed  and  accuracy  of  computational  quantum  chemistry  (QC)  have  allowed  for 
the  calculation  of  molecular  geometries  and  energies  for  small  molecules  with  relative  ease. 
Nevertheless,  it  is  still  impractical,  due  to  computational  expense,  to  simulate  a  large  system  of  reacting 
molecules  (100  or  more)  over  a  relevant  time  scale  using  a  purely  quantum  mechanical  description.  To 
simulate  large  systems  at  an  atomistic  scale,  it  is  typical  to  use  a  force  field  approach  instead  of  QC.  For 
example,  the  MM3  force  field  [15]  has  been  used  successfully  in  the  description  of  hydrocarbons,  and 
more  general  force  fields,  such  as  DREIDING  [16],  are  able  to  provide  relatively  accurate  results  for  a 
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number  of  different  elements.  However,  to  allow  for  the  simulation  of  large  systems,  these  force  fields 
use  sufficiently  simple  descriptions  of  bonding  and  are  not  able  simulate  chemical  reaction  because  they 
cannot  model  the  bond  breaking  and  formation.  One  approach  to  analyzing  such  systems  is  through  the 
development  of  reactive  force  fields. 

Recently,  a  reactive  ReaxFF  force  field  has  been  developed  to  study  the  oxidation  of  ammonia 
borane  (NH3BH3)  [17].  ReaxFF  is  a  reactive  force  field  method  which  uses  “bond-order”  to  describe 
bonds  as  a  function  of  inter-atomic  distance,  as  opposed  to  explicitly  defined  bonding  [18].  This  allows 
for  continuous  breaking  and  forming  of  bonds,  as  the  bond  order  is  updated  at  every  iteration  as  atom 
distance  changes.  Furthermore,  a  ReaxFF  force  field  can  be  parameterized  using  data  derived  from  either 
experimental  or  QC  calculations,  so  that  no  extensive  experimental  data  is  needed  [18].  The  ReaxFF 
method  has  been  discussed  in  depth  in  previous  publications,  which  should  be  referred  to  for  greater 
detail  [18-20].  . 

Using  the  ReaxFF  reactive  force  field  developed  previously  [17]  with  data  acquired  through 
density  functional  theory  (DFT)  calculations,  this  paper  discusses  the  development  of  a  chemical  kinetic 
mechanism  based  on  MD  simulation  results.  The  first  step  was  to  conduct  MD  simulations  and  analyze 
the  results  to  identify  chemical  pathways.  From  these,  the  intermediate  species  and  elementary  reactions 
of  importance  can  be  deduced.  For  species  for  which  there  were  no  thermochemical  data  available,  the 
enthalpy  and  entropy  of  formation  were  calculated  using  QC  methods,  along  with  specific  heat  as  a 
function  of  temperature.  Next,  the  rates  for  new  reactions  were  estimated.  In  this  work,  simple  collision 
theory  was  used  to  attain  initial  estimates  of  reaction  rates.  Chemical  kinetics  calculations  were  then 
performed  and  compared  with  the  MD  results.  Refinements  to  the  mechanism  are  made  by  performing 
additional  MD  simulations  to  isolate  specific  steps  in  the  reaction  progress  and  by  using  more  accurate  ab 
initio  reaction  rate  calculation  methods  for  the  most  important  reactions,  for  example,  as  determined  by 
sensitivity  analysis  on  the  chemical  kinetics  calculations. 

The  MD  simulations  used  to  study  the  chemical  pathways  through  which  AB  oxidation  proceeds 
were  all  performed  under  NVT  conditions  (constant  number  of  atoms,  volume,  and  temperature).  A 
Berendsen  thermostat  [21]  with  a  0.1  ps  damping  constant  was  used  to  control  temperature.  For  an 
accurate  simulation,  the  force  field  must  be  able  to  model  flexible  molecules  with  flexible  bonds  that 
exhibit  translational,  rotational,  torsional,  and  vibrational  motion.  Normally,  a  time  step  one  order  of 
magnitude  smaller  than  the  shortest  motion,  or  approximately  0.5  to  1.0  fs,  is  required.  A  smaller  time 
step  is  required  when  using  ReaxFF  because  the  charges  and  bond  orders  are  allowed  to  change  at  every 
time  step.  For  ReaxFF  MD  simulations  which  exceed  temperatures  of  2500  K,  a  time  step  of  0.1  fs  is 
preferred  for  efficient  coverage  of  the  phase  space  and  collisions  and  reactions  to  occur  smoothly. 

For  this  work,  molecules  were  simulated  inside  of  a  2.5  nm  periodic  cube.  The  use  of  high 
pressure  and  temperature  is  required  in  ReaxFF  simulations  in  order  to  create  enough  successful 
molecular  collisions  within  reasonable  simulation  times.  For  the  simulations  discussed  here,  the  initial 
system  pressure  in  approximately  170  bar,  however  at  the  highest  temperatures  it  can  exceed  2500  bar. 
The  requirement  of  high  pressure  makes  MD  simulations  best  suited  for  modeling  systems  such  as  rocket 
and  gun  propellants  or  shock  induced  reactions. 

To  model  the  chemistry  of  a  system,  the  thermochemical  parameters  of  each  species  must  be 
known,  including  heat  of  formation,  entropy  of  formation,  and  specific  heat.  Since  borane  amines  have 
been  studied  recently  for  chemical  vapor  deposition  of  BN  films  and  hydrogen  storage,  thermochemical 
data  have  been  calculated  for  many  of  the  pertinent  B/N/H  species  with  ab  intio  QC  methods,  as  reported 
in  the  literature  [22-24].  However,  the  thermochemical  properties  of  B/N/H/O  species  are  not  available. 
Therefore,  DFT  calculations  were  undertaken  to  calculate  these  parameters  for  this  study.  The  DFT- 
calculations  used  to  derive  the  thermochemical  properties  were  performed  using  the  B3LYP  functional 
[25].  This  method  combines  exact  HF  exchange  with  the  Becke  generalized  exchanged  functional  of 
Becke  [26]  and  Lee,  Yang,  and  Parr  [27]  correlation  functional.  For  the  B,  H,  O,  and  N  atoms,  all 
electrons  were  included  using  a  modified  variant  of  the  Pople  6-311G**++  basis  set  [28,  29].  All 
calculations  were  performed  using  the  Jaguar  7.8  program  package  [30]. 
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The  first  step  in  calculating  the  heat  of  formation  at  the  standard  temperature  (AH^sk)  for  a 
given  molecule  (M)  was  to  calculate  the  atomization  energy  ED0,  as  given  by  Eq. 

z  (1) 

atoms 


where  e0(X)  is  the  total  energy  of  atom  X,  v  is  the  stoichiometric  coefficient  for  a  given  atom,  s0(M)  is  the 
total  energy  of  molecule  M,  and  szpe(M)  is  the  zero  point  energy  of  molecule  M. 

Given  the  atomization  energy  of  the  molecule,  along  with  the  heats  of  formation  for  the 
constituent  atoms,  the  heat  of  formation  of  molecule  M  at  0  K  was  computed  according  to  Eq.  . 
Subsequently,  the  heat  of  formation  at  298  K  was  computed  using  Eq. 

A Hf  (M,0K)=Yj  wA Hf  (X)  -  Z£>0  (2) 

atoms 

A Hf  (M, 29SK)  =  A Hf  ( M,0K )  +  [Hm (298 K)  -  Hu (0/Q] 

-Y,o[Hx(298K)-Hx(0K)\  (3) 

atoms 


Using  vibrational  analysis,  the  specific  heat  and  entropy  of  the  molecule  was  calculated  for  a 
specified  temperature.  By  generating  data  for  a  range  of  temperatures,  the  parameters  for  a  fourth  order 
curve  fit  of  specific  heat  were  calculated  for  two  temperature  ranges  (200-1000  K  and  1000-6000  K),  as 
well  as  integration  constants  describing  enthalpy  and  entropy,  to  be  used  in  equilibrium  and  kinetics 
programs  such  as  CHEMKIN  [31]. 

Equilibrium  calculations  were  performed  using  CHEETAH  4.0  [32]  with  the  JCZS  product 
library  developed  by  Hobbs  and  Baer  [33].  Chemical  kinetic  calculations  were  performed  using 
CHEMKIN,  release  4.0  [31].  To  describe  hydrogen  oxidation  and  ammonia  chemistry,  the  95  reaction 
kinetic  mechanism  published  by  Lindstedt  et  al  [34]  was  used.  The  B/O/H  reactions  were  described  using 
the  first  30  reactions  in  the  mechanism  published  by  Yetter  et  al [13].  For  some  of  the  important 
reactions  in  AB  decomposition,  reaction  rates  have  been  proposed  by  Nguyen  et  al  [35]  using  ab  initio 
QC  calculations  and  RRKM  theory.  At  1  atm,  the  Arrhenius  reaction  rate  equations  for  AB  hydrogen 
elimination  and  B-N  bond  breaking  are  given  by  Eq.  and  ,  respectively. 

h3nbh3  ->  h2nbh2  +  h2 

£  =  3.841£  +  12-exp(-35-6%r|  (4) 

H3NBH3^NH3+BH3 

A;  =  1.545£  +  12-exp(-21-9%r)  (5) 

Nguyen  et  al  also  studied  the  catalytic  effect  of  BH3  on  hydrogen  elimination  from  AB.  QC  calculations 
show  that  a  collision  with  BH3  can  significantly  reduce  the  energy  barrier  for  hydrogen  release. 

Therefore,  a  separate  rate  equation  for  this  reaction  was  specified  [35],  as  given  by  Eq.  . 

H3NBH3+BH3  — » H2NBH2  +H2+  BH3 

A:  =  6.330£  +  llexp(-22-4%r)  (6) 
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For  a  preliminary  evaluation  of  reaction  rate  parameters,  simple  collision  theory  provides  a 
simple,  upper-limit  estimation.  Although  this  method  of  estimation  is  subject  to  a  larger  degree  of 
uncertainty,  it  accomplishes  the  objective  of  establishing  an  initial  set  of  rate  constants  which  are 
physically  reasonable  and  collectively  self-consistent.  According  to  simple  collision  theory,  the 
frequency  factor  A  may  be  written  as  in  Eq.,  where  p  is  the  steric  factor,  which  is  the  ratio  of  successful 
collisions  to  total  collisions,  oAB  is  the  hard  sphere  collision  diameter  for  a  collision  between  molecule  A 
and  molecule  B,  kb  is  Boltzmann’s  constant,  T  is  temperature,  and  //^is  the  reduced  mass. 


/  \  1/2 

%nkbT' 

Bab  y 


BaB 


mAmB 
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° AB  ~ 
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The  collision  diameter  where  crA  and  crB  are  the  hard  sphere  diameters  of  the  colliding  molecules.  Hard 

sphere  diameters  for  well  known  species  were  taken  from  Kee  et  al  [36],  while  those  not  available  in  the 
literature  (most  boron  containing  species)  were  estimated  using  a  semi-empirical  relationships  first 
proposed  by  Hirschfelder  and  Eliason  [37].  The  general  form  for  this  relationship  is  given  by  Eq.  ,  where 
T  is  the  mean  atom  radius  and  A  and  B  are  constants.  Hirschfelder  and  Eliason  [37]  proposed  a 
relationship  for  diatomics  where  A=4/3  and  B=1.8  A.  Svehla  [38]  extended  the  relationship  to  triatomics 
with  constants  of  A=10/9  and  B=1.7  A.  Based  on  values  for  a  in  the  literature  [36,  38],  this  formula  was 
extended  to  larger  molecules  (four  to  eight  atoms)  for  this  study,  where  constants  of  A=8/23  and  B=2.95 
A  were  found  to  match  well  with  data  for  a  wide  range  of  larger  molecules.  Values  for  T  were  taken 
from  Clementi  et  al  [39]. 


'N 
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One  particularly  important  reaction  not  previously  studied  is  the  second  hydrogen  elimination 
step  from  AB.  To  estimate  the  reaction  rate  for  this  step,  a  more  accurate  approach  was  taken.  Using  a 
transition  state  search  in  Jaguar  7.8  with  B3LYP/6-311G**++,  a  transition  state  for  hydrogen  elimination 
from  H2NBH2  was  found,  as  pictured  in  Figure  1.  The  reaction  rate  was  then  calculated,  using  the  initial 
and  transition  state  description,  in  the  ChemRate  program  [34].  An  Arrhenius  rate  equation  was 
determined,  as  shown  in  Eq.. 

^ h2nbh2^>hnbh+h2  ~  5.488£T  +  13exp|  ^3 ^/rt)  (^) 


Figure  1:  Reaction  coordinate  of  hydrogen  e 


imination  from  H2NBH2,  showing  the  reactant, 


transition  state,  and  product,  as  calculated  with  Jaguar  7.8  using  B3LYP/6-311G**++ 


Figure  2  shows  the  equilibrium  products  as  a  function  of  temperature  for  a  stoichiometric  reaction 
of  AB  and  molecular  oxygen  at  a  pressure  of  1  atm.  Over  the  range  of  temperatures  shown,  the  main 
equilibrium  products  always  consist  of  gaseous  H20,  N2,  and  to  a  lesser  extent,  H2.  However,  the  species 
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containing  boron  atoms  changes  significantly  over  this  range.  Below  650  K,  the  majority  of  boron  exists 
as  solid  phase  HOBO.  Between  650  and  900  K,  boron  exists  mostly  as  gaseous  B3H3O6  combined  with 
gaseous  B3H303  and  liquid  phase  B203.  Liquid  B203  is  the  primary  boron  containing  species  from  900  K, 
until  gaseous  HOBO  overtakes  it  around  1700  K. 


Figure  2:  Mole  fractions  of  equilibrium  products  of  a  stoichiometric  reaction  of  AB  and  molecular 

oxygen  at  a  pressure  of  1  atm 

Figure  3  (a.)-(c.)  shows  snapshots  from  a  temperature  ramping  NVT  MD  simulation  starting  with 
20  AB  and  45  02  molecules  in  a  2.5  nm  periodic  cube,  with  an  imposed  heating  rate  of  0.00522  K/fs.  The 
composition  was  determined  as  a  stoichiometic  mixture  according  to  Eq.  .  It  should  be  noted  that  the 
stoichiometric  mixture  ratio  remains  the  same  regardless  of  whether  the  boron  forms  B203  or  HOBO. 

NH3BH3+^02  +HOBO  +  ^H20  (10) 
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Figure  3:  Snapshots  from  a  temperature  ramp  simulation  with  20  AB  molecules  and  45  02 
molecules,  (for  color  figure,  green=boron,  blue=nitrogen,  white=hydrogen,  red=oxygen)  showing 
a.)initial  minimized  geometry  b.)T=2500  K  c.)  T=5000  K 

Figure  3  (a.)  shows  the  initial  geometry  with  a  2.5  nm  transparent  periodic  cube,  where  the  molecules 
have  been  randomly  placed  and  subsequently  energy  minimized.  Figure  3  (b.)  shows  a  snapshot  at  the 
middle  of  the  simulation,  where  the  thermostat  is  set  to  2500  K,  and  Figure  3  (c.)  shows  the  final  frame  of 
the  simulation  with  thermostat  setting  of  5000  K.  At  the  point  in  the  simulation  which  Figure  3  (b.) 
pictures,  the  AB  molecules  have  all  been  dehydrogenated  leaving  no  AB  molecules  and  a  significant 
amount  of  molecular  hydrogen  in  the  system.  Most  of  this  H2  has  not  started  to  react  with  the  oxygen,  but 
there  are  some  H20  and  H02  molecules  in  the  system.  Several  of  the  dehydrogenated  AB  molecules  have 
been  attacked  by  oxygen,  forming  HNBOH  molecule.  In  Figure  3  (c.),  it  can  be  seen  that  many  of  the 
expected  equilibrium  products  have  been  formed.  Molecular  hydrogen  has  been  converted  to  H20,  the 
majority  of  nitrogen  exists  as  N2,  and  the  boron  exists  primarily  as  HOBO. 

Figure  4  shows  a  plot  of  the  major  species  observed  versus  temperature  for  the  temperature 
ramping  NVT  MD  simulation  pictured  in  Figure  3.  Consistent  with  MD  simulations  on  AB  in  the 
absence  of  oxygen,  some  of  the  AB  is  converted  to  NH3  and  BH3  through  B-N  bond  cleavage  in  the  low 
temperature  range.  The  NH3  and  BH3  do  not  oxidize,  and  the  B-N  bonds  begin  to  reform  as  the 
temperature  increases  beyond  1250  K.  Around  1000  K,  H2  begins  to  be  generated  through  the  uni- 
molecular  hydrogen  elimination  from  AB.  Just  as  shown  previously  in  the  snapshots,  H2  molecules  are 
oxidized  to  form  H20,  while  the  B-N  containing  molecules  undergo  a  reaction  sequence,  which  ultimately 
ends  up  producing  N2  and  HOBO. 
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Figure  4:  Temperature  ramping  MD  simulation  of  AB  and  oxygen  molecules 

Figure  5  is  a  flow  diagram  describing  the  reaction  sequence  through  which  AB  progresses  to  form 
the  equilibrium  products  H20,  N2,  HOBO,  and  possibly  B203  in  the  MD  simulations.  As  mentioned 
previously,  the  B-N  bond  in  AB  is  relatively  weak  and  tends  to  break  at  low  temperatures,  but  as  their 
concentrations  become  greater  the  reverse  reaction  quickly  becomes  dominant  and  AB  is  reformed.  The 
reaction  sequence  of  the  AB  molecule  is  down  the  left  side  of  the  diagram,  undergoing  a  first  and  then 
second  hydrogen  elimination,  contributing  to  a  large  concentration  of  H2  in  the  system,  and  resulting  in 
the  HNBH  molecule.  The  strain  energy  for  the  HNBH  molecule  to  undergo  hydrogen  elimination  is  too 
great,  and  the  B-N  bond  in  this  molecule  is  quite  strong  (-170  kcal,  as  determined  by  a  coordinate  scan  in 
Jaguar  7.8  using  B3LYP/6-311G**++).  For  the  reaction  to  proceed  further  towards  the  equilibrium 
products,  the  boron  side  of  this  molecule  will  react  with  02  or  radicals  formed  from  H2/02  chemistry, 
since  the  boron  atom  is  electron  deficient  and  the  B-H  bond  is  the  weakest  in  the  molecule.  The  MD 
simulation  results  show  the  HNBOH  molecule  to  be  a  major  product.  The  B-N  bond  is  still  relatively 
strong  in  the  HNBOH  molecule,  with  a  bond  dissociation  energy  of  143  kcal/mol.  However,  the  B-N 
bond  must  be  severed  in  the  conversion  to  the  equilibrium  products.  Investigating  the  minor  species 
formed  in  the  MD  simulations  reveals  that  the  boron  side  of  HNBOH  is  still  likely  to  react  with  OH,  O, 
and  02,  while  the  nitrogen  side  lowers  its  energy  by  forming  another  bond  with  a  hydrogen,  leading  to  the 
formation  of  what  turns  out  to  be  a  very  low  energy  species,  H2NB(OH)2.  The  B-N  bond  energy  for 
H2NB(OH)2  has  been  reduced  to  approximately  111  kcal,  making  it  the  weakest  bond  in  the  molecule. 
The  unfavorable  entropy  of  this  molecule  means  that  it  easily  dissociates  at  high  temperatures.  B-N  bond 
cleavage  results  in  the  formation  of  an  NH2  an  B(OH)2.  The  B(OH)2  molecule  can  easily  form  HOBO 
with  a  hydrogen  extraction  by  the  NH2.  The  pathways  for  N/H  species  to  N2  are  well  studied  and  can  be 
described  by  reactions  found  in  the  literature  [34].  Furthermore,  depending  on  the  system  temperature,  it 
may  be  favorable  for  HOBO  to  be  converted  to  B203  and  H20,  the  pathways  of  which  have  also  been 
previously  reported  [13]. 
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Figure  5:  Flow  chart  describing  the  pathways  of  AB  oxidation,  formulated  through  ReaxFF  MD 

simulations  of  AB  and  O2. 


Figure  6:  Molecular  geometries  of  B/N/H/O  species,  calculated  using  Jaguar  7.8  [24]  with  the 

B3LYP  functional  and  6-311G**++  basis  set 

Figure  6  contains  the  molecular  geometries  of  B/N/H/O  species  that  were  deemed  relevant  to  the 
kinetics  of  AB  oxidation  from  the  MD  simulations  and  for  which  prior  data  was  not  available.  Table  1 
reports  the  atomization  energy,  heat  of  formation  and  entropy  of  formation  for  molecules  shown  in  Figure 
6,  as  well  as  previously  studied  molecules  relevant  to  this  study,  as  calculated  with  DFT  calculations 
using  the  B3LYP  functional  and  the  6-31 1G**++  basis  set.  Although  it  is  known  that  B3LYP  is  not  the 
most  accurate  method  for  calculating  energies  (post  hartree-fock  methods  such  as  CCSD  are  generally 
considered  to  have  higher  accuracy  for  energy  calculations),  it  is  a  very  reliable  and  efficient  method  and 
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gives  satisfactory  results  for  this  initial  study.  As  mentioned  above,  the  H2NB(OH)2  molecule  exhibited  a 
very  low  heat  of  formation  (-179.90  kcal/mol  at  298  K),  the  lowest  of  the  molecules  examined  here. 
However,  equilibrium  calculations  still  show  this  molecule  is  not  a  favored  product  due  to  the  effect  of 
entropy. 

Table  1:  Thermochemical  parameters  calculated  using  Jaguar  7.8  [30]  with  the  B3LYP  functional 
and  6-311G**++  basis  set  and  Eq.  compared  to  values  from  the  literature  "[41]  b  [23] c  [22] 


calculated  literature 


Molecule 

Molecular  Energy 
[Hartrees] 

Zero  Point  Energy 
[kcal] 

ID0 

[kcal] 

AHf(OK) 

[kcal/mol] 

AHf(298K) 

[kcal/mol] 

S(298K) 

[cal/mol*K] 

AHf(298K) 

[kcal/mol] 

S(298K) 

[cal/mol*K] 

HOBO 

-176.0402181 

12.75 

427.59 

-122.98 

-123.65 

58.04 

-133.99a 

57.3a 

nh3 

-56.576047 

21.52 

273.55 

-6.56 

-8.24 

48.16 

-10.973 

46.07a 

bh3 

-26.621103 

16.40 

267.91 

22.56 

25.49 

48.55 

25.50a 

44.90a 

B203(gas) 

-275.620437 

11.87 

637.89 

-191.04 

-191.26 

69.47 

-199.533 

67.82a 

H3NBH3  (gas) 

-83.25008 

43.70 

568.90 

-11.45 

-15.85 

59.43 

-13.5b,-19.95c 

58. 3C 

H2NBH2 

-82.07418 

30.04 

475.01 

-21.09 

-23.84 

55.91 

-15.9b,-23.2c 

54.78c 

HNBH 

-80.83477 

14.89 

342.75 

7.64 

12.05 

37.27 

12.8b,  11.0C 

41.62c 

HNBOH 

-156.13090 

18.90 

469.67 

-60.44 

-61.43 

62.41 

- 

- 

H2NBO 

-156.166813 

20.06 

491.05 

-81.82 

-83.10 

61.28 

- 

- 

H2NB(0H)2 

-232.696329 

37.86 

747.79 

-176.19 

-179.90 

69.13 

- 

- 

HNBO 

-155.477472 

11.55 

382.17 

-24.70 

-25.28 

59.80 

- 

- 

H2NBOH 

-156.708865 

27.11 

508.97 

-47.98 

-50.21 

63.59 

- 

- 

HNB(OH)2 

-232.002435 

29.12 

636.27 

-116.45 

-119.28 

69.47 

- 

- 

Figure  7  shows  the  results  of  a  constant  pressure,  constant  energy  chemical  kinetics  calculation 
performed  in  CHEMKIN  4.0  [31].  The  initial  mixture  consisted  of  a  stoichiometric  mixture  of  AB  and 
air.  The  mechanism  contains  the  reactions  from  the  literature  as  discussed  earlier  [13,  34],  as  well  as  52 
new  reversible  reactions,  which  are  listed  in  Table  SI  in  supplemental  material.  At  very  early  times,  not 
visible  in  Figure  7,  the  B-N  bond  breaks  in  AB,  quickly  inundating  the  system  with  NH3  and  BH3.  As  in 
the  MD  simulations,  the  reverse  reaction  quickly  becomes  important  and  the  NH3  and  BH3  are  converted 
back  to  AB.  At  the  same  time,  the  AB  hydrogen  elimination  reaction  takes  place,  converting  all  of  the 
AB  to  H2NBH2  and  H2  after  approximately  1  ms.  The  H2  is  slowly  consumed  through  reaction  with  02, 
injecting  some  radicals  into  the  system  and  forming  H20,  which  leads  to  an  increased  temperature.  After 
approximately  20  ms,  H2NBH2  begins  to  be  converted  to  HNBH  and  more  H20.  The  mole  fraction  of 
HNBH  steadily  increases  until  just  before  the  rapid  oxidation  event  at  62  ms,  where  it  is  converted  to 
HOBO  and  N2.  A  sensitivity  analysis  shows  that  the  most  important  elementary  reaction  in  the  overall 
reaction  progress  is  the  chain  branching  H  +  02^>  OH  4-  O  ,  while  reactions  that  compete  for  the 

consumption  of  H  atoms,  such  as  HBNH2  +  H  — »  HNBH  +  H2  and  H2BNH  +  H  — »  HNBH  +  H2 , 
are  the  most  inhibiting  to  the  ignition  of  the  system. 
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Figure  7:  Plot  of  species  mole  fraction  and  temperature  versus  time  for  a  constant  pressure  and 
energy  kinetics  calculation  performed  in  Chemkin  4.0 

Due  to  its  high  hydrogen  content,  AB  has  the  potential  to  increase  performance  if  used  as  a  fuel  in 
chemical  propulsion  applications.  Molecular  dynamics  simulations  have  been  performed  to  examine  the 
oxidation  of  AB  using  a  previously  developed  ReaxFF  reactive  force  field  developed  to  describe  B/N/H/O 
chemistry  [11].  The  results  of  these  simulations  were  used  to  identify  chemical  pathways,  important 
intermediate  species,  and  a  series  of  elementary  reactions  for  the  conversion  of  AB  and  molecular 
hydrogen  to  equilibrium  species.  DFT  calculations  were  performed  using  the  B3LYP  functional  and  the 
6-311G**++  basis  set  in  the  Jaguar  7.8  software  to  calculate  thermochemical  properties  of  unknown 
species.  Estimated  reaction  rates  were  combined  with  known  reactions  to  assemble  a  chemical  kinetic 
mechanism  for  AB  oxidation.  This  mechanism  was  used  to  perform  constant  pressure  and  enthalpy 
kinetics  calculations.  These  calculations  show  the  reaction  proceeds  as  follows:  rapid  decomposition  of 
AB  (H3NBH3)  to  form  hydrogen  and  H2NBH2,  reaction  of  the  hydrogen  to  form  radicals  and  H20  and 
raise  the  temperature,  dehydrogenation  of  H2NH2  to  form  HNBH,  and  finally,  an  ignition  event  which 
rapidly  converts  the  intermediates  to  the  equilibrium  products  . 

These  results  will  be  used  to  design  experiments  to  further  investigate  and  validate  the  reaction 
kinetics  of  AB  oxidation.  Through  refinement  procedures,  critical  reactions  in  this  mechanism  will  be 
identified  and  targeted  for  more  accurate  reaction  rate  estimate  using  ab  initio  methonds,  such  as 
transition  state  theory. 

B.  Confined  Rapid  Thermolysis  Studies  of  Ammonia  Borane  (with  Stef  Thynell,  Penn 

State  University) 

Contemporary  AB  research,  recently  reviewed  by  Staubitz  et  al.  [42],  is  primarily  concerned  with 
generating  hydrogen  for  use  in  proton-exchange  fuel  cells,  and  using  catalysis  or  other  methods  to  provide 
increased  hydrogen  release  at  low  temperatures.  The  thermal  decomposition  of  AB,  first  studied  by  Hu  et 
al.  [43],  has  recently  been  investigated  in  a  series  of  publications  emanating  from  a  group  at  Friedberg 
University  of  Mining  and  Technology  [44,  45,  46].  These  studies  rely  mostly  on  combination  of 
TGA/DSC  experiments. 

The  first  of  these  studies  [44]  shows  that  for  sufficiently  slow  heating  rates  (0.05  K/min),  or 
isotherms  above  70°C  for  sufficient  time  (-30  hours),  AB  fully  completes  the  first  decomposition  step 
below  its  melting  temperature  (104°C)  [44],  demonstrating  that  at  the  low  limit  of  heating  rate,  AB  can 
fully  undergo  its  first  exothermic  decomposition  step  without  changing  phase.  Following  the  lengthy 
isotherms  ranging  from  70-90  °C,  the  resulting  condensed-phase  product  was  determined  to  be  polymeric 
aminoborane  [(NH2BH2)X]  [44].  At  increased  heating  rates,  the  compound  will  melt  prior  to  the  bulk  of 
the  decomposition,  with  gas  release  resulting  in  a  vigorous  foaming. 

Additionally,  the  investigations  of  Wolf  et  al  [44]  and  Baitalow  et  al.  [45]  at  slow  heating  rates 
(<  5  K/min)  show  that  the  percent  mass  loss  increased  with  an  increasing  heating  rate,  while  the  volume 
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of  evolved  gas  (as  measured  by  water  displacement)  was  independent  of  heating  rate.  Assuming 
hydrogen  was  the  only  product  released  into  the  gas  phase,  Baitalow  et  al  measured  approximately  2.2 
mol  H2/mol  of  AB  when  heating  up  to  190°C,  corresponding  to  a  total  mass  loss  of  14.3  percent  for  all 
heating  rates  [45].  At  a  heating  rate  of  1  K/min,  the  gravimetric  mass  loss  recorded  was  roughly  15 
percent,  corresponding  closely  with  their  volumetric  measurement  of  hydrogen.  However,  at  a  heating 
rate  of  5  K/min,  the  gravimetric  mass  loss  was  approximately  23  percent,  about  1.6  times  the  mass  of  the 
measured  hydrogen  [45].  This  result  from  Baitalow  et  al  [45]  demonstrates  that  a  certain  percentage  of 
the  decomposition  gases  are  either  condensable  at  room  temperature  or  are  water-soluble.  Furthermore, 
the  amount  of  these  condensable/water-soluble  products  increases  with  increasing  heating  rate.  Gas- 
phase  FTIR/MS  studies  performed  in  conjunction  with  these  TGA/DSC  studies  give  evidence  that  these 
products  consist  of  monomeric  aminoborane  (BH2NH2),  borazine  (c-[BHNH]3),  and  small  amounts  of 
other  BNHX  species  [45]. 

Stowe  et  al  [47]  have  studied  the  mechanism  of  thermal  decomposition  of  AB  below  its  melting 
point  by  probing  the  solid  phase  using  nuclear  magnetic  resonance  (NMR)  spectroscopy.  Using  the  NMR 
data,  they  have  proposed  a  mechanism  through  which  the  AB  can  react  to  form  polymeric  amino-borane 
(PAB)  and  release  hydrogen  into  the  gas  phase.  The  NMR  spectra  show  the  reaction  is  initiated  with  the 
formation  of  a  “mobile”  phase  of  AB,  which  then  leads  to  the  formation  of  an  ionic  compound,  the 
diamoniate  of  diborane  or  DADB  [(NH3)2BH2]+[BH4]‘.  It  is  suggested  that  the  more  reactive  DADB  will 
react  with  the  surrounding  AB,  resulting  in  hydrogen  release  and  amino-borane  polymerization. 

Another  detailed  study  of  AB  decomposition  has  recently  been  published  by  Frueh  et  al  [13],  in 
which  materials  characterization  techniques  are  used  in  addition  to  TGA/DSC  analysis  to  understand  the 
chemical  pathways  of  AB  decomposition.  This  study  claims  to  observe  a  third  hydrogen  release  at  a 
temperature  of  1170°C,  with  significantly  slower  reaction  rate  than  the  first  two  hydrogen  releases.  The 
resulting  material  after  prolonged  heating  was  determined  to  be  hexagonal  boron  nitride  through  the  use 
of  several  material  characterization  techniques  [48]. 

For  energetic  materials,  the  kinetic  decomposition  pathways  may  be  highly  dependent  on  the 
heating  rate.  For  example,  Oyumi  and  Brill  found  that  with  RDX  decomposition,  changing  the  external 
heating  rate  results  in  a  shift  from  CH20  and  N20  as  the  main  decomposition  products  at  low  heating 
rates  to  HCN  and  N02  at  high  heating  rates  [49].  Therefore,  to  understand  the  decomposition  behavior  of 
propellant  ingredients,  it  is  widely  recognized  that  the  process  must  be  studied  using  both  slow  and  rapid 
heating  conditions.  Previous  works  [43-46]  have  only  examined  decomposition  at  slow  heating  rates 
(<10  K/min),  where  the  time  scale  for  reaction  is  on  the  same  order  as  the  time  scale  for  heating.  In  this 
work,  the  decomposition  of  AB  is  studied  using  confined  rapid  thermolysis  (CRT).  In  this  experiment, 
the  sample  is  heated  at  rates  of  approximately  2000  K/sec  to  a  set  temperature  where  diagnostics  of  the 
isothermal  decomposition  are  performed.  TGA/DSC  experiments  are  also  performed  to  investigate 
changes  in  mass  loss  and  heat  release  with  increased  heating  rate. 

The  ammonia  borane  used  in  this  study  was  purchased  from  Boroscience  International.  The 
material  was  specified  to  be  98%  pure  by  the  supplier  using  proton  and  nB  NMR  spectroscopy.  Its 
physical  appearance  is  that  of  a  white  powder. 

The  CRT  experiment  was  previously  described  in  [50].  In  this  experiment,  a  small  amount  of 
material  (-0.5  mg)  is  heated  in  a  narrow  space  between  two  flat  surfaces.  Two  copper  rods,  53  mm  in 
height  and  15.6  mm  in  diameter,  are  heated  to  a  preset  isothermal  condition,  using  high-watt  density 
cartridge  heaters.  The  copper  rods  are  then  compressed  against  the  sample  holder,  consisting  of  a 
retaining  ring  and  thin  metal  foils,  using  a  pneumatic  piston  cylinder.  The  sample  holder  maintains  a 
sealed  300  pm  gap  between  the  heaters,  except  for  a  rectangular  slit  that  allows  decomposition  gases  to 
escape  into  the  cool  atmosphere  and  into  the  modulated  beam  of  the  FTIR.  The  CRT  apparatus  is  used  in 
conjunction  with  two  diagnostic  tools  to  examine  the  gas-phase  species:  rapid  scanning  FTIR 
transmission  spectroscopy,  and  time-of-flight  mass  spectrometry  (ToF-MS).  These  two  diagnostic  tools 
complement  one  another  and  are  most  effective  when  used  in  conjunction. 

For  the  CRT/gas-phase  FTIR  spectroscopy,  the  setup  is  contained  in  a  constant  pressure  chamber 
filled  with  N2.  The  modulated  FTIR  beam,  which  enters  and  exits  the  chamber  via  ZnSe  windows,  passes 
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through  the  decomposition  gases  that  emanate  from  a  rapidly  heated  sample.  The  Bruker  IFS  66/S  FTIR 
spectrometer  is  programmed  to  acquire  the  transmissivity  spectrum  with  a  spectral  resolution  of  2  cm'1 
and  a  temporal  resolution  of  50  ms. 

As  previously  described  by  Chowdhury  et  al.  [51],  the  ToF-MS  measurements  take  place  with  an 
identical  CRT  heating  system  inside  a  low-pressure  chamber.  A  commercially  available  ToF-MS 
spectrometer  is  used.  Several  stages  of  vacuum  are  used  to  generate  a  molecular  beam  from  the  evolved 
decomposition  gases.  The  gases  first  pass  through  a  100  pm  orifice  into  the  first  stage,  pumped  down  to 
10'4  Torr.  Subsequently,  the  gases  pass  through  a  1  mm  Ni  skimmer  to  the  second  stage,  pumped  to  10"6 
Torr,  and  a  0.5  mm  slit  to  a  third  stage,  pumped  to  10'7  Torr.  Next,  an  electron  beam,  with  an  impact 
ionization  energy  set  at  70  eV,  is  incident  on  the  molecular  beam  for  a  period  of  3  ps.  The  resulting 
positive  ions  created  by  the  interaction  of  the  two  beams  are  accelerated  across  a  known  potential  (250  V 
between  grid  plates,  2.5  kV  total  potential).  Since  the  potential  and  distance  traveled  are  well  known,  the 
mass-to-charge  ratio  is  found  by  measuring  the  “time  of  flight.” 

AB  decomposition  results  in  a  significant  amount  of  condensed-phase  products.  To  identify  and 
characterize  condensed-phase  species,  an  FTIR/KBr  pellet  technique  was  used.  A  small  amount  (-0.5 
mg)  of  AB  was  heated  at  a  known  temperature  for  a  known  period  of  time  in  the  CRT  apparatus  described 
above.  The  resulting  condensed-phase  products  remaining  on  the  metal  foil  were  then  mixed  with  -150 
mg  of  KBr  powder.  The  resulting  mixture  was  then  pressed  into  a  pellet.  The  pellet  is  placed  inside  a 
housing,  purged  with  N2,  and  aligned  with  the  modulated  FTIR  beam.  The  same  FTIR  spectrometer 
acquires  the  transmissivity  spectrum,  where  the  spectrum  for  a  pure  KBr  pellet  serves  as  the  background. 
In  addition  to  examining  the  residue  left  on  the  metal  foil  after  heating,  a  glass  slide  was  used  to  collect 
products  in  the  gaseous  product  stream  which  can  easily  condense  out  at  room  temperature.  A  diagram 
showing  the  position  of  the  glass  slide  is  pictured  in  8.  After  being  collected,  the  condensate  is  removed 
from  the  glass  slide  by  grinding  with  KBr  powder.  The  sample  is  then  pressed  into  a  KBr  pellet,  for  later 
FTIR  analysis. 


Cold  glass  slide 
collects  gaseous 
species  which  may 
condense  out 

Rectangular  slit  allows 
gaseous  decomposition 
species  to  escape 


Figure  8  Diagram  showing  the  position  of  glass  slide  for  collecting  condensate  from  gaseous 

decomposition  products 


Thermal  analysis  was  conducted  using  a  simultaneous  TGA/DSC  apparatus  (Netzch  STA-449 
Jupiter),  which  has  a  silicon  carbide  furnace  with  a  maximum  temperature  of  1550°C  and  a  maximum 
heating  rate  of  50  K/min.  The  samples  were  placed  in  alumina  crucibles  covered  with  thin  alumina  lids 
that  had  a  small  hole  to  release  the  evolved  gases.  Sample  masses  were  approximately  2.5  mg  for  each 
test.  The  sample  mass  was  limited  by  the  amount  of  material  swelling  during  decomposition. 

The  sample  chamber,  which  is  vacuum-tight,  was  evacuated  and  backfilled  with  an  argon 
atmosphere  prior  to  each  measurement,  and  was  continuously  purged  with  argon  during  the  experiments. 
A  top-loading,  low-drift  balance  with  a  resolution  of  1  pg  measures  the  change  in  sample  mass.  Two  S- 
type  thermocouples  are  attached  to  a  platinum  sample  carrier,  which  holds  both  a  reference  and  sample 
crucible,  to  measure  heat  flow  from  the  sample.  The  onset  melt  temperatures  of  pure  metal  samples  of 
indium,  tin,  zinc,  aluminum  and  gold  were  used  to  calibrate  the  temperature  measurements.  Prior  to  each 
experiment,  a  baseline  measurement  is  made  to  compensate  for  buoyancy  and  specific  heat  effects  in  the 
TGA  and  DSC  measurements,  respectively. 

Figure  9  shows  a  plot  of  signal  intensity  for  different  mass-to-charge  (m/z)  ratios,  as  measured 
by  CRT/ToF-MS.  Also  included  in  this  figure  are  the  species  assignments  corresponding  to  each  mass- 
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to-charge  ratio.  The  first  species  to  appear  in  the  signal  are  H2+  (m/z=2)  and  NH3+  (m/z=17),  followed 
shortly  after  by  m/z=28  attributed  to  either  NH210BH2+,  NH2nBH+,  or  NHUBH2+.  Subseqently,  strong 
peaks  are  shown  at  m/z=27  and  29,  which  may  correspond  to  NH10BH2  or  NHnBH+  and  NH2UBH2+, 
respectively.  Mass-to-charge  ratios  of  27,  28,  and  29  can  all  be  attributed  to  the  fragmentation  spectrum 
of  molecular  H2NBH2.  At  significantly  later  times,  signals  at  m/z=53,  80,  and  81  are  detected,  suggesting 
the  presence  of  borazine  (c-[NHBH]3).  At  a  thermolysis  temperature  of  180°C,  the  AB  is  consumed  very 
quickly  indicated  by  the  early  decrease  in  the  evolution  of  ammonia. 

— • —  m/z=2  — • —  m/z=27  — * —  m/z=29  • —  m/z=80 


— ■ —  m/z=17  — * —  m/z=28  — * —  m/z=53  — *—  m/z=81 


Figure  9:  Mass  to  charge  ratio  signals  for  ToF-MS  of  decomposition  products  from  AB  at  180°C 

Figure  10  shows  a  series  of  FTIR  spectra  from  the  CRT  experiment  at  160°C  with  0.6  mg  of  AB. 
The  first  spectrum,  taken  at  650  ms,  shows  two  very  strong  peaks  at  965  and  930  cm'1  and  some  slightly 
weaker  peaks  at  1083  and  993  cm'1,  indicating  the  presence  of  ammonia  (NH3).  In  the  second  spectrum, 
taken  at  750  ms,  the  ammonia  peaks  have  increased  in  strength,  and  quite  a  few  new  peaks  have  also 
appeared.  Peaks  in  the  region  of  3300-3500  cm'1  are  indicative  of  N-H  stretching  frequencies  and  the 
peaks  at  3479  and  3451  cm'1  are  attributed  to  the  asymmetric  and  symmetric  N-H  stretches  of  H2NBH2. 
The  peaks  from  2400-2600  cm'1  can  be  attributed  to  B-H  stretching.  The  peaks  at  2597  and  2493  cm"1 
match  well  with  the  asymmetric  and  symmetric  B-H  stretching  frequencies  reported  in  the  literature  [52, 
53].  Multiple  peaks  for  each  B-H  stretching  mode  likely  arise  due  to  the  presence  of  two  stable  boron 
isotopes  (nB  and  10B)  in  relatively  high  quantities.  Likewise,  the  frequencies  for  B-N  stretching  can  be 
attributed  to  peaks  at  1368  cm'1  (10B-N)  and  1340  cm"1  (nB-N),  as  they  match  well  with  data  from  the 
literature  for  the  H2NBH2  spectrum  [52].  The  peak  at  1625  cm"1  also  matches  well  with  values  reported  in 
the  literature  [52,  53],  and  it  is  attributed  to  -NH2  scissoring. 
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Figure  10:  FTIR  transmission  spectra  of  the  gaseous  decomposition  products  from  a  0.6  mg  sample 
of  AB  heated  to  160°C  in  the  CRT  experiment. 
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Figure  1 1 :  FTIR  transmission  spectra  of  the  residue  from  AB  pyrolysis  at  160°C  subsequently 
heated  to  250°C  in  the  CRT  experiment 

The  third  spectmm,  taken  at  900  ms,  shows  the  ammonia  peaks  declining  while  the 
aforementioned  H2NBH2  peaks  stay  at  about  the  same  intensity.  New  peaks  are  apparent  in  this  spectmm 
at  wavenumbers  of  3493,  2555,  and  1465  cm"1,  indicating  the  presence  of  a  new  substance.  The  intensity 
of  the  H2NBH2  peaks  decreases  in  the  subsequent  spectra,  and  at  1750  ms,  they  have  nearly  vanished. 
The  prominent  peaks  in  the  final  spectra  (3492,  2520,  1465,  and  917  cm"1)  all  correspond  with  the  very 
strong  peaks  observed  by  Niedenzu  et  al.  [54]  for  the  gas-phase  IR  spectmm  of  borazine  (c-(NH-BH)3). 

Figure  11  shows  the  FTIR  transmission  spectra  of  the  gas-phase  products  resulting  from  the 
residue  of  the  previously  described  experiment  (0.6  mg  sample  of  AB  rapidly  heated  to  160°C)  being 
rapidly  heated  to  250°C  in  the  CRT  experiment.  The  first  spectmm  shown,  taken  at  500  ms,  has  peaks 
which  provide  evidence  that  both  aminoborane  (H2NBH2)  and  borazine  (c-(NH-BH)3)  are  produced.  The 
peaks  at  1368  and  1340  cm"1,  also  evident  in  the  previous  spectra,  are  attributed  to  the  10B-N  and  nB-N 
stretches  in  H2NBH2.  The  very  strong  peak  at  1465  cm"1  is  attributed  to  B-N  stretching  in  the  borazine 
molecule.  The  feature  at  2520  cm"1  is  a  combination  of  the  N-H  stretching  peaks  for  both  molecules, 
while  the  peaks  at  3480  and  3493  cm"1  are  attributed  to  the  B-H  stretching  for  H2NBH2  and  borazine, 
respectively.  A  prominent  peak  attributed  to  -NH2  symmetric  bending  of  H2NBH2  is  observed  at  1618 
cm"1.  The  other  strong  peaks  in  the  spectmm,  at  1597  and  917  cm"1,  were  previously  observed  in  the 
FTIR  spectmm  of  borazine  [54].  The  second  spectmm  in  the  sequence,  taken  100  ms  later,  shows  all  the 
peaks  in  the  spectmm  growing  in  intensity.  No  significant  new  peaks  are  observed.  Subsequent  spectra 
taken  at  1000  and  1250  ms  show  the  peaks  for  H2NBH2  diminishing,  followed  by  a  decrease  in  the 
borazine  peaks  in  the  final  spectra.  It  is  of  note  that  no  peaks  corresponding  to  ammonia  appear  in  the 
spectra,  suggesting  that  AB  is  not  present  in  the  residue. 

Figure  12  shows  a  sequence  of  FTIR  transmission  spectra  taken  of  the  condensed-phase  residue 
which  remained  on  the  foil  after  AB  thermolysis  in  the  CRT  experiment.  The  spectmm  labeled  “0  sec” 
corresponds  to  unheated  AB  as  received  from  Boroscience  International.  The  subsequent  spectra  were 
acquired  from  the  condensed-phase  residue  of  decomposed  AB  after  being  heated  for  a  known  period  of 
time.  For  all  of  the  spectra  in  ,  the  sample  was  heated  at  140°C. 

The  absorption  bands  in  the  spectra  from  3000-3500  cm"1  can  be  attributed  to  N-H  stretching 
vibration.  An  expanded  view  of  this  region  is  shown  in  Figure  13.  The  peaks  in  the  neat  AB  spectmm  at 
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3253  and  3194  cm'1  correspond  well  with  previously  reported  wavenumbers  for  asymmetric  and 
symmetric  N-H  stretching  for  AB  in  a  KBr  pellet  [55].  After  heating  at  140°C  for  4  sec,  the  original 
peaks  appear  to  still  be  present,  with  the  symmetric  N-H  stretch  slightly  shifted.  However,  a  new  peak 
has  appeared  at  3305  cm'1.  After  further  heating,  the  peaks  from  the  neat  AB  diminish  compared  to  the 
new  feature,  resulting  in  a  broad  absorption  curve  with  a  peak  at  3291  cm'1  after  12  sec.  Previous  studies 
have  observed  absorption  peaks  attributed  to  N-H  stretching  in  polymers  amino-borane  in  this  range 
(3300  cm'1  by  Komm  et  al  [56]  for  (H2NBH2)X,  3301  cm'1  by  Boddeker  et  al  [57]  for  (H2NBH2)5).  For 
amine  boranes,  the  frequency  of  N-H  vibration  increases  with  increasing  B-N  bond  order,  as  Frueh  et  al 
[48]  have  noted.  This  provides  evidence  that  the  AB  has  dehydrogenated  and  formed  PAB. 


Figure  12:  Overlayed  time  sequence  of  FTIR  spectra  from  condensed  phase  AB  residue  from  CRT 
experiments  at  140°C.  Residue  was  mixed  with  KBr  and  pressed  into  a  pellet. 


Wavenumbers  (cm-1) 


Figure  13:  N-H  stretching  frequencies  in  AB  and  AB  residue  for  various  heating  times  at  140°C 

Figure  14  shows  the  absorption  bands  attributed  to  B-H  stretching  between  2000  and  2500  cm'1. 
For  the  neat  AB,  asymmetric  and  symmetric  B-H  stretching  frequencies  are  observed  at  2348  and  2288 
cm'1,  slightly  higher  than  those  reported  in  the  literature  (2330  and  2280  cm"1  [55]).  After  just  4  sec  of 
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heating,  the  spectral  features  have  changed  significantly  in  this  region.  New  peaks  are  observed  at  2385, 
2331,  2290,  2236,  and  2222  cm"1.  These  peaks  continue  to  make  up  the  spectra  after  8  and  12  sec  of 
heating.  These  peaks  seem  to  correspond  best  with  those  reported  by  Baumann  et  al.  [46]  for  (BNHX)Z 
polymer  (decomposition  products  of  (H2NBH2)X),  suggesting  the  AB  has  undergone  more  than  one  H2 
elimination. 


Figure  14:  B-H  stretching  frequencies  in  AB  and  AB  residue  for  various  heating  times  at  140°C 
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Figure  15:  Fingerprint  region  of  the  FTIR  transmission  spectrum  for  AB  and  AB  residue  from 
samples  heated  at  140°C  for  varying  amounts  of  time 

Figure  15  shows  an  enlarged  view  of  the  FTIR  spectrum  from  1650-1000  cm"1,  in  which  the 
upper  part  of  what  is  often  referred  to  as  the  “fingerprint”  region.  In  this  region,  there  is  usually  present 
the  fundamental  vibrational  modes  due  to  molecular  bending.  The  neat  AB  shows  peaks  at  1608  cm"1, 
and  1384  cm"1,  which  correspond  well  with  values  reported  in  the  literature  for  NH3  deformation  [55,  58]. 
The  peak  at  1422  cm"1  in  the  neat  AB  has  not  been  reported  in  the  literature,  but  is  tentatively  attributed  to 
N-H  rocking.  The  peaks  at  1170  and  1162  cm"1  closely  match  wavenumbers  in  the  literature  for  BH3 
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deformation  [55,  58],  as  does  the  peak  at  1066  for  the  NBH  rock  in  AB.  Spectra  of  the  AB  residue  all 
show  a  peak  at  1566  cm'1,  which  closely  matches  the  spectral  location  reported  by  Boddeker  et  al.  for 
NH2  deformation  in  the  cyclic  (H2NBH2)2  [57].  The  peak  at  1422  cm'1  vanishes  when  the  AB  is  heated, 
and  the  peak  at  1384  cm'1  shifts  to  the  left  and  is  located  at  approximately  1389  cm'1.  It  is  probable  that 
this  peak  is  no  longer  related  to  NH3  deformation,  but  instead  can  be  attributed  to  B-N  stretching,  which 
occurs  near  1400  cm'1  for  double  B-N  bonds  [46,  48,  57].  The  features  in  neat  AB  at  1 170  and  1162  cm'1 
are  broadened  after  heating  and  result  in  a  peak  shifted  to  1165  cm'1  after  12  sec.  A  peak  at  1126  cm'1 
appears  upon  heating  and  grows  relative  to  the  other  peaks  at  longer  heating  time.  This  peak  could  be 
related  to  BH2  bending,  and  a  similar  peak  was  reported  by  Boddeker  et  al.  [57]  for  cyclic  (H2BNH2)4. 
Lastly,  it  appears  that  the  peak  attributed  to  NBH  rocking  is  still  present  after  heating,  and  after  12  sec  it 
has  shifted  slightly  to  1064  cm'1.  It  should  be  noted  that  some  small  shifts  in  band  location  may  be 
attributed  to  changes  in  the  association  of  molecules  in  the  condensed  phase. 

Neat  AB 


Figure  16:  FTIR  transmission  spectrum  showing  B-H  stretching  frequencies  measured  from  neat 
AB,  condensate  collected  on  a  glass  slide,  and  residue  removed  from  sample  foils  after  16  seconds  at 


120°C 


NeatAB 


Figure  17:  FTIR  transmission  spectrum  showing  fingerprint  region  for  neat  AB,  condensate 
collected  on  a  glass  slide,  and  residue  removed  from  sample  foils  after  16  seconds  at  120° 
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Figure  16  shows  the  FTIR  spectrum  in  the  B-H  stretching  region  for  the  condensate  collected  on 
the  glass  slide,  as  compared  with  neat  AB  and  the  residue  from  rapid  thermolysis  at  120°C  for  16 
seconds.  The  peaks  in  the  condensate  correspond  closely  with  those  in  the  residue,  especially  those  at 
2226,  2292,  and  2442  cm'1.  The  shoulder  which  appears  in  the  residue  around  2380  cm'1  is  notably 
absent  in  the  condensate,  and  the  strongest  absorption  peak  is  shifted  to  higher  wavenumbers,  closer  to 
that  found  in  the  neat  AB.  Figure  17  shows  the  fingerprint  region  for  the  condensate,  neat  AB,  and 
thermolysis  residue  from  heating  at  120°C  for  16  seconds.  The  spectrum  corresponding  to  the  condensate 
again  matched  well  with  the  condensed  phase  residue  of  AB  after  heating  at  120°C  in  this  region.  As  in 
the  thermolysis  residue  ,  the  spectrum  of  the  condensate  showed  the  shifting  of  the  N-H  rocking 
absorption  peak  at  1422  cm'1  to  1408  cm'1.  Both  of  these  spectra  also  showed  the  appearance  of  peaks  at 
1568  cm"1,  1188  cm"1,  and  1100  cm"1.  Many  of  the  peaks  from  neat  AB  are  still  observed  in  both  the 
residue  and  condensate.  These  results  suggest  that  the  condensate  has  roughly  the  same  composition  as 
the  condensed  phase  residue  at  early  times.  They  both  appear  to  be  comprised  of  some  AB  along  with 
dehydrogenated  amino-borane  polymers.  The  implications  of  this  result  are  that  H2NBH2  in  the  gaseous 
products  will  readily  condense  out  on  a  cold  surface  and  begin  to  polymerize.  Also,  this  result  suggests 
that  some  AB  will  sublimate  into  the  gas  phase,  but  it  decomposes  too  fast  to  be  observed,  even  by  the 
fast  diagnostics  used  in  the  gas  phase  CRT  experiments. 
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Figure  18:  Thermo-gravimetric  and  differential  scanning  calorimetry  data  for  AB  heated  at  a  rate 
of  20  K/min  in  an  Argon  atmosphere 

Figure  1 8  shows  data  typical  of  a  TGA/DSC  experiment  with  ammonia  borane.  The  case  plotted 
in  is  from  a  test  with  a  heating  rate  of  20  K/min.  Two  distinct  mass  losses  appear  in  Figure  18.  The  first, 
with  an  onset  temperature  of  approximately  1 13°C,  corresponds  to  a  mass  loss  of  9.4  percent  of  the  initial 
mass.  The  second  mass  loss  corresponds  to  19.0  percent  of  the  initial  mass,  and  has  an  onset  temperature 
of  roughly  152°C.  AB  samples  were  heated  up  to  1400°C;  however,  the  signal  output  from  either  the 
TGA  or  DSC  instrument  do  not  reveal  any  noteworthy  changes  above  200°C.  This  is  inconsistent  with 
the  work  of  Frueh  et  al  [48],  who  report  a  slow  mass  loss  starting  at  1 170°C.  The  shape  of  the  mass  loss 
curve  corresponds  well  with  those  reported  by  Baitalow  et  al  [45]  for  lower  heating  rates  (1  and  5 
K/min).  Also,  in  agreement  with  the  work  of  Baitalow  et  al ,  the  percent  mass  loss  increases  with 
increased  heating  rate,  is  a  scatter  plot  of  the  percent  mass  loss  as  a  function  of  heating  rate  with  data 
from  this  study  and  Baitalow  et  al  [45].  These  data  show  that  mass  loss  is  nearly  logarithmic  function  of 
heating  rate  over  the  range  of  2-50  K/min. 
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To  build  on  previous  studies  of  AB  decomposition  at  low  heating  rates  [44,  45,  46,  48],  we  have 
investigated  decomposition  at  high  heating  rates.  Knowledge  of  decomposition  pathways  at  high  heating 
rates  is  important  for  studying  the  combustion  of  propellant  ingredients.  The  TGA  mass  loss  data  shows 
that  the  mass  loss  continues  to  increase  logarithmically  with  heating  rates  up  to  50  K/min.  The  increased 
mass  loss  is  presumable  due  to  a  greater  amount  of  condensable  or  water-soluble  decomposition  products 
as  observed  by  Baitalow  et  al  [45]  for  lower  heating  rates. 

CRT/FTIR/ToF-MS  experiments,  where  samples  are  heated  at  rates  of  approximately  2000  K/sec 
to  an  isothermal  condition,  were  performed  to  further  investigate  the  effect  of  high  heating  rates  on  the 
decomposition  products.  The  data  from  the  ToF-MS  reveal  that  di-hydrogen  (H2),  ammonia  (NH3), 
monomeric  aminoborane  (H2BNH2),  and  borazine  (c-[NHBH]3)  emerge  into  the  gas  phase.  The  FTIR 
transmission  spectroscopy  measurements  support  the  presence  of  these  species.  Aside  from  ammonia, 
Baitalow  et  al.  [45]  identified  all  these  species  at  low  heating  rates.  Ammonia  is  likely  produced  from  the 
cleavage  of  the  B-N  bond  in  AB,  as  shown  in  Eq.  .  This  is  a  coordinate  bond,  which  has  a  relatively  weak 
dissociation  energy  (31.3  kcal/mol  [59]).  The  cleavage  of  this  bond  has  been  shown  to  be  important  at 
low  temperatures  in  MD  simulations  [60].  After  careful  examination  of  the  acquired  FTIR  spectra, 
borane  (BH3)  has  not  been  detected  in  the  gas-phase.  This  is  likely  due  to  its  high  reactivity,  suggesting 
that  it  takes  part  in  condensed-phase  reactions  rather  than  escaping  into  the  gas  phase.  By  contrast,  the 
fairly  un-reactive  NH3  much  more  readily  escapes  the  solid  or  liquid  matrix. 

NH3BH3  ->  NH3  +  bh3  (1  r) 

The  production  of  gaseous  ammonia  at  high  heating  rates  may  be  part  of  the  reason  for  increase 
mass  loss  at  high  heating  rates.  A  larger  effect  is  likely  related  to  the  rate  of  formation  of  high  molecular 
weight  polymers.  As  noted  by  Wolf  et  al  [44],  the  condensed-phase  products  following  the  first 
hydrogen  release  are  mostly  polymeric  aminoborane  or  PAB  ([NH2BH2]X).  With  low  temperature 
isotherms  (>90°C)  or  low  heating  rates,  the  AB  molecule  releases  an  H2  molecule  prior  to  its  melting 
temperature,  and  the  resulting  NH2BH2  molecule  is  able  to  polymerize  with  like  molecules.  Stowe  et  al 
[47]  have  suggested  this  process  is  initiated  by  the  formation  of  the  ionic  compound  DADB.  The 
activation  energy  for  such  polymerization  processes  is  typically  negative.  Thus,  for  a  slow  heating  rate, 
where  temperature  will  remain  lower  for  longer  times,  polymerization  is  enhanced.  At  higher  heating 
rates,  AB  will  melt  prior  to  decomposing.  Furthermore,  since  heating  is  much  more  rapid  than  the 
reaction,  the  decomposition  will  occur  at  a  higher  temperature,  where  polymerization  is  slower  and  less 
favorable  thermodynamically.  As  the  time  scale  for  polymerization  increases,  more  monomeric 
aminoborane  will  be  able  to  diffuse  into  the  gas  phase,  therefore  increasing  the  amount  of  mass  lost  to  the 
gas  phase  at  higher  heating  rates  and  higher  temperatures.  The  collection  of  condensate  with  the  glass 
slide  shows  that  gaseous  monomeric  aminoborane  will  condense  out  and  polymerize  on  a  cold  surface. 
The  spectra  of  the  condensate  also  show  evidence  that  some  AB  may  be  present,  even  though  it  is  not 
detected  with  the  fast  gas  phase  diagnostics.  This  suggests  that  some  AB  may  sublime,  but  it  is  rapidly 
converted  to  H2NBH2  and  H2  in  the  gas  phase. 
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Figure  19:  Percent  mass  loss  versus  heating  rate  for  TGA  experiments  with  AB 

Hydrogen  elimination  from  PAB  has  been  reported  to  occur  above  130°C  [46],  and  condensed 
phase  polyimino  borane  or  PIB  ([HNBH]X)  is  known  to  be  formed  as  a  result.  Furthermore,  the  gas  phase 
FTIR  measurements  in  this  study  show  that  monomeric  aminoborane  and  borazine  evolved  from  the 
residue  of  AB  decomposition,  likely  consisting  of  PAB  and  PIB.  The  pathways  for  the  conversion  of  PAB 
and  PIB  to  borazine  have  been  studied  by  Nutt  and  McKee  [61].  Figure  20  shows  a  schematic  of  the 
pathways  for  AB  decomposition.  As  indicated,  the  formation  of  PAB  is  favored  at  low  heating  rates, 
while  higher  heating  rates  favor  gaseous  monomeric  aminoborane. 

HNBH  +  H 

2  2  2 


Figure  20:  Schematic  of  condensed  phase  thermolysis  of  AB  and  evolved  gas  phase  products. 

Since  propellant  ingredients  will  experience  rapid  heating  in  practical  applications,  it  is  important 
to  study  the  decomposition  pathways  under  similar  heating  rates,  as  these  processes  can  be  rate 
dependent.  Confined  rapid  thermolysis  was  used  to  investigate  the  decomposition  of  AB  at  high  heating 
rates.  In  this  experiment,  AB  samples  were  rapidly  heated  to  a  preset  temperature  level,  where 
decomposition  was  examined  under  isothermal  conditions.  This  type  experiment  differs  from  that  of 
TGA/DSC,  where  the  heating  rate  is  slow  and  the  temperature  varies  continuously.  Studies  of  the  gas- 
phase  products  by  FTIR  spectroscopy  and  by  ToF-MS  show  the  evolution  of  H2  and  NH3  followed  a  short 
time  after  by  H2NBH2.  At  later  times,  borazine  (c-N3B3H6)  is  observed.  Aside  from  NH3,  these  products 
have  all  been  previously  observed  to  evolve  from  decomposition  at  low  heating  rates  [45].  While  the 
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amount  of  hydrogen  evolved  has  been  shown  to  remain  the  same,  the  total  mass  of  evolved  products 
increases  when  the  sample  is  more  quickly  heated  to  a  temperature  above  its  decomposition  point. 
Condensate  collected  on  a  cold  glass  slide  in  the  stream  of  gaseous  decomposition  products  has  an  FTIR 
spectrum  with  matches  closely  with  that  of  polymeric  aminoborane  ([H2NBH2]X). 


III.  Modeling  and  Simulation  of  Thermo-Chemical  Behavior  of  Nano-Energetic 
Materials  (Georgia  Tech) 

A.  Molecular  Dynamics  Simulations  of  Physicochemical  Behavior  of  Aluminum, 

Nickel,  and  Nickel-Coated  Nanoaluminum  Particles 

Aluminum  particles  are  extensively  used  in  many  propulsion  and  energy-conversion  applications, 
due  to  their  favorable  energetic  properties  [62].  They  are  covered  by  an  amorphous  oxide  layer,  which  is 
typically  0.5-4  nm  thick  [63]. The  specific  thickness  of  the  oxide  layer  depends  on  the  temperature  of  the 
particle  and  the  duration  of  exposure  to  the  oxidizing  environment.  When  the  particle  is  exposed  to  the 
oxidizing  gas  for  a  sufficient  period,  the  oxide  layer  thickness  saturates  at  a  value  of  4  nm  [64].  The  oxide 
layer  protects  the  particle  from  further  attack  by  the  oxidizer  molecules.  Aluminum  particles  with 
diameter  greater  than  100  pm  ignite  only  upon  melting  of  the  oxide  layer  at  2350  K  [65].  This  results  in  a 
long  ignition  delay  and  a  slow  rate  of  energy  release.  It  is  thus  desirable  to  reduce  the  ignition  temperature 
of  micron-sized  aluminum  particles. 

The  formation  of  the  oxide  layer  on  aluminum  particles  can  be  significantly  inhibited  by  applying 
transition  metal  coatings,  which  have  higher  melting  points  than  the  aluminum  [5].  For  example,  the  bulk 
melting  point  of  nickel  is  1728  K.  Foley  et  al.  [66]  synthesized  nano-aluminum  powders  by  thermal 
decomposition  of  an  alane-adduct  solution  in  the  presence  of  a  titanium  catalyst  under  an  inert 
atmosphere.  The  resulting  material  was  used  to  reduce  complexes  of  gold,  nickel,  palladium,  and  silver. 
The  reduction  process  yielded  nano-aluminum  particles  coated  with  transition  metals.  Typically,  an  oxide 
layer  of  the  transition  metal  is  formed,  which  can  have  thickness  on  the  order  of  few  nm  [66-67].  Among 
all  transition  metals  considered  in  Foley  et  al.’s  study,  nickel  resulted  in  the  highest  active  aluminum 
content.  Encapsulation  of  aluminum  particles  with  a  nickel  shell  results  in  lower  ignition  temperatures, 
due  to  inter-metallic  reactions  between  aluminum  and  nickel  atoms  [68].  For  example,  the  ignition 
temperature  of  a  2.38  mm  aluminum  particle  in  air  decreases  from  2350  to  13 13  K  when  the  oxide  layer  is 
replaced  with  a  nickel  coating.  Temperature  runaway  is  also  observed  in  inert  environments  [68].  The 
substitution  of  a  nickel  coating  for  the  oxide  coating  increases  the  flame  speed  of  an  aluminum  particle 
dust  in  air  by  a  factor  of  1.5-4  [69].  Nickel-coated  aluminum  particles  can  potentially  be  employed  to 
manufacture  nickel  aluminides  (NixAli_x),  which  are  attractive  for  use  in  aircraft  turbines  due  to  their  good 
strength  characteristics,  excellent  corrosion  and  oxidation  resistance,  and  high  melting  points  [70].  A 
better  understanding  of  its  thermo-chemical  behavior  is,  thus,  helpful  for  propulsion  and  material 
synthesis  applications. 

Nano-sized  particles  have  unusual  energetic  properties  due  to  their  high  percentage  of  surface 
atoms  and  the  excess  energy  associated  with  these  atoms  [71-73].  As  the  particle  size  decreases  from  30 
to  3  nm,  the  percentage  of  surface  atoms  increases  from  5  to  50  %  [74].  The  melting  temperature  of  a 
nano-aluminum  particle  increases  from  473  K  at  2  nm  to  a  bulk  value  of  937  K  at  approximately  8  nm 
[75].  The  ignition  temperature  can  be  as  low  as  933  K,  possibly  due  to  polymorphic  phase 
transformations  in  the  oxide  layer  and/or  melting  of  the  aluminum  core  [63,76].  The  oxide  layer, 
however,  occupies  a  significant  portion  of  the  particle  mass  at  nano-scales.  For  example,  a  38  nm 
aluminum  particle  has  -47.5  wt.  %  oxide  [77].  For  combustion  applications,  therefore,  it  is  desirable  to 
replace  the  oxide  coating  with  such  favorable  metallic  coatings  as  nickel.  Nickel  atoms  can  participate  in 
oxidation  and  intermetallic  reactions,  thereby  resulting  in  a  significant  reduction  in  the  ignition  delay. 
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Molecular  dynamics  (MD)  simulations  can  be  used  to  provide  insight  into  the  melting  and 
ignition  characteristics  of  nano-scale  materials.  Puri  and  Yang  [75]  performed  MD  simulations  to  analyze 
the  effect  of  particle  size  on  melting  of  nascent  aluminum  particles.  The  thermo-mechanical  behavior  of 
oxide  coated  aluminum  particles  has  also  been  analyzed  [78].  Molecular  dynamics  simulations  on  nickel- 
coated  aluminum  particles  have  thus  far,  however,  been  limited.  Most  efforts  were  made  to  study  the 
possibility  of  fragmentation  of  the  shell  caused  by  core  melting  at  a  fixed  core/particle  size  [79,  80],  with 
little  attention  paid  to  diffusion  and  inter-metallic  reactions  [80,  81]. 

Systematic  studies  of  the  effect  of  particle  size  on  the  physiochemical  processes  after  melting  of 
the  core  are  yet  to  be  performed.  In  the  present  study,  the  thermochemical  behavior  of  nickel-coated 
aluminum  particles  is  investigated  via  MD  simulations  over  a  wide  range  of  temperatures  to  characterize 
melting,  diffusion,  and  inter-metallic  reactions.  The  core  diameter  of  the  particles  of  interest  is  in  the 
range  of  3-12  nm,  and  shell  thicknesses  of  0.5,  1,2  and  3  nm  are  considered.  Special  attention  is  paid  to 
quantify  the  effect  of  the  core  diameter  and  nickel  shell  on  relevant  physiochemical  phenomena. 

In  the  present  study,  both  isobaric-isothermal  (NPT)  and  isochoric-isoenergetic  (NVE)  ensembles 
are  employed.  NPT  ensemble  is  used  to  study  the  thermochemical  behavior  of  the  particle  in  the  presence 
of  external  heating  at  a  constant  pressure  condition,  while  NVE  ensemble  is  used  to  investigate  self¬ 
heating  of  the  particle  due  to  inter-metallic  reactions  under  adiabatic  conditions.  In  NPT  ensemble,  a 
system  of  N  atoms  is  coupled  to  an  external  source  by  introducing  additional  variables  into  the 
Lagrangian.  Assuming  that  atoms  behave  as  classical  point-like  masses,  the  Lagrangian,  Z,  is  expressed  as 
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where  Mis  a  constant  fictitious  mass  associated  with  the  volume  of  the  system,  m  the  mass  of  the  atom,  qt 
the  generalized  coordinate,  V  the  volume,  treated  as  a  dynamic  variable,  P  the  pressure,  U  the  potential 
energy,  T  the  temperature,  Q  the  inertia  factor,  g  a  parameter,  kB  the  Boltzmann  constant,  and  5  the  degree 
of  freedom  of  the  thermostat.  In  this  method,  volume  is  treated  as  a  dynamic  variable  [82]  and  the 
parameter  5  controls  heat  exchange  between  the  system  and  the  reservoir  [83].  Substituting  the 
Lagrangian  in  the  Euler-Lagrange  equation  [23],  the  equations  of  motion  take  the  form 
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where  qt  and  q.  denote  the  first-order  and  second-order  derivatives  of  qt  with  respect  to  time,  and  Ft  is 

the  net  force  on  atom  i.  In  NVE  ensemble,  the  Lagrangian  is  the  difference  between  the  kinetic  and 
potential  energies  of  the  set  of  atoms.  The  resulting  equation  of  motion  is  given  by 

dU 

miqi=-—.  (14) 

dq, 

The  above  system  of  equations  is  numerically  integrated  using  a  fifth-order  predictor-corrector 
algorithm.  The  time  step  is  chosen  as  1  fs  because  the  time  scale  of  vibration  of  atoms  is  of  the  same 
order.  Negligible  improvement  in  the  model  results  is  obtained  when  the  time  step  is  reduced  from  1  to 
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0.1  fs.  The  heating  rate  is  another  important  adjustable  parameter.  Shibuta  and  Suzuki  [85]  analyzed  the 
effect  of  the  cooling  rate  on  the  solidification  of  metal  nano-particles  using  MD  simulations.  The 
solidification  temperature  decreased  with  increasing  cooling  rate.  In  the  present  work,  a  parametric  study 
was  conducted  to  determine  the  optimum  value  of  the  heating  rate  in  the  range  10'3-10_1  K/fs.  The  result  is 
shown  in  Fig.  21.  The  time  step  places  restriction  on  the  minimum  value  of  the  heating  rate.  The  analysis 
reveals  that  a  heating  rate  lower  than  10"2  K/fs  increases  the  total  computational  time  dramatically  with 
only  little  change  in  the  model  results,  while  a  higher  value  leads  to  substantially  different  results.  As  a 
result,  a  heating  rate  of  10"2  K/fs  is  adopted  in  the  present  study. 


Figure  21:  Variation  of  potential  energy  of  7  nm  (10976  atoms)  aluminum  particle  with 

temperature  at  different  heating  rates 


To  calculate  macroscopic  properties  of  the  system  from  the  positions  and  momenta  of  all  atoms,  an 
ensemble  average  is  required  [84].  The  pressure  is  calculated  using  the  virial  equation  of  state,  as  a 
function  of  temperature  and  forces  experienced  by  all  the  atoms  [86]. 


P  =  NkBT  + 


(15) 


where  rt  is  the  position  vector  of  atom  i.  The  melting  point  of  the  system  is  identified  based  on  the 
variations  in  the  potential  energy,  Lindemann  index,  thermal  displacement,  and  translational-order 
parameter.  The  Lindemann  index,  S ,  is  a  measure  of  the  root-mean-square  fluctuations  of  the  inter-atomic 
distance,  and  is  given  by  [87] 


S  = 


(16) 


It  is  expected  to  increase  abruptly  during  melting.  The  translational-order  parameter,  X,  is  expressed  as 

[88] 
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where  k  is  an  arbitrary  vector  of  the  lattice  and  Z  is  the  side  length  of  the  simulation  domain.  In  solids, 
atoms  simply  translate  about  their  lattice  positions,  and  hence,  there  is  a  perfect  order  in  the  crystal.  Such 
order,  however,  disappears  in  a  liquid  state,  and  the  corresponding  translational-order  parameter  is 
reduced  by  an  order  of  magnitude.  The  thermal  displacement  is  defined  as 


Su  = 


(18) 


where  the  subscripts  i  and  0  refer  to  atom  i  and  the  initial  state,  respectively.  Upon  melting,  the 
aluminum  atoms  diffuse  into  the  nickel  shell.  The  onset  of  diffusion  is  marked  by  the  sudden 
increase  in  the  core  radius,  rc  [79] 


Z=1 


(19) 


where  rcm  is  the  position  vector  of  the  center  of  mass. 

Potential  Function 

To  close  the  formulation,  the  potential  energy  of  the  system  of  atoms  needs  to  be  specified.  For 
metals,  the  majority  of  the  potential  functions  are  based  on  the  embedded  atom  method  [89].  In  this 
method,  the  energy  needed  to  embed  an  atom  in  the  electron  gas  is  included  to  capture  the  physics  of 
metallic  bonding.  One  such  potential  function,  which  is  formulated  by  Cleri  and  Rosato  [90],  has  been 
fitted  to  the  structural  and  thermodynamic  properties  of  bulk  aluminum  and  nickel.  It  takes  the  following 
form: 
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where  V  is  the  pair-potential  function  and  p  the  electron  density  function.  The  first  term  describes  the 
electrostatic  interactions  between  the  atoms,  while  the  second  term  accounts  for  the  energy  needed  to 
embed  an  atom  in  the  electron  gas.  The  pair-potential  and  electron  density  functions  are  expressed  as 
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Table  2  Parameters  for  potential  function  [91,  92] 


Ni-Ni 

Al-Al 

Ni-Al 

A(eV) 

0.0376 

0.1221 

0.0597 

C(eV) 

1.0700 

1.3160 

1.2898 

P 

16.9990 

8.6120 

15.714 

q 

1.1890 

2.5160 

1.1550 

r<>  (A) 

2.4910 

2.8637 

2.5001 

The  parameters  in  the  potential  function  are  given  in  Table  2  [91,  92].  Previous  studies  indicate  that  the 
potential  function  can  predict  the  structure,  phase  stability,  and  thermodynamic  properties  of  Ni-Al  alloys 
(especially,  B2-NiAl  and  Ni3Al)  with  reasonable  accuracy.  Furthermore,  considerable  success  has  been 
achieved  in  simulations  thermal  behaviors  of  aluminum  particles  coated  with  Ni  and  Ni3Al  using  this 
potential  function.  As  a  result,  it  is  employed  in  the  present  study. 

Pure  Aluminum  and  Nickel 

The  theoretical  framework  is  employed  to  calculate  the  structural  and  thermodynamic  properties  of 
aluminum  and  nickel.  An  FCC  lattice  is  adopted  to  calculate  the  initial  position  vectors  for  a  known 
number  of  atoms.  For  bulk  materials,  a  periodic  boundary  condition  is  enforced  in  all  the  three  spatial 
directions.  A  free-surface  boundary  condition  is,  however,  prescribed  for  nano-particles.  The  crystal  is 
allowed  to  equilibrate  at  300  K  prior  to  the  heating  simulation. 

Bulk  Materials 

The  properties  of  surface-free  bulk  systems  of  aluminum  and  nickel  are  first  treated.  The  lattice 
constant  and  cohesive  energy  are  calculated  by  equilibrating  the  crystal  in  an  NVE  ensemble.  To 
determine  the  melting  point  and  latent  heat  of  melting,  the  crystal  is  heated  externally  in  an  NPT 
ensemble.  Table  3  shows  calculated  values  of  the  cohesive  energy,  lattice  constant,  and  latent  heat  of 
melting  and  compares  them  with  experimental  data.  The  discrepancy  is  less  than  2  %,  demonstrating  the 
accuracy  of  the  potential  function.  Figure  22  shows  the  melting  of  4000-atom  bulk  aluminum  in  vacuum. 
Melting  is  observed  at  1060  K,  at  which  point  sharp  variations  occur  in  the  translational-order  parameter 
and  Lindemann  index.  Such  a  trend  is  characteristic  of  the  structural  melting  phenomenon,  which  is 
observed  for  surface-free  materials  [96].  Note  that  the  fluctuations  that  are  typically  observed  in  MD 
simulations  decrease 

Table  3  Cohesive  energy,  lattice  constant,  and  latent  heat  of  melting  of  bulk  materials 


Aluminum 


Nickel 


Model  Experiment93"95  Model  Experiment93'95 


Ec  (eV/atom)  -3.337  -3.340  -4.437  -4.435 

a  (A)  4.048  4.050  3.491  3.520 

AHm  (kJ/mol)  10,610 _ 10.470 _ 17.360 _ 17.160 
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Figure  22  Translational  order  parameter  (X)  and  Lindemann  index  (8)  as  a  function  of  temperature 
showing  the  melting  of  bulk  aluminum  and  nickel  in  vacuum 


with  increasing  number  of  atoms.  The  system  is  superheated  above  the  experimental  melting  point  of  bulk 
material  by  127  K,  which  can  be  attributed  to  the  absence  of  a  nucleation  site  for  melting.  The  structural 
melting  point  is  typically  greater  than  the  thermodynamic  counterpart  by  20  %  [96,97].  The  resulting 
asymptotic  value  of  the  melting  temperature  of  aluminum  particles  is  883  K,  which  differs  from  the 
experimental  value  by  about  5  %.  A  similar  pattern  is  observed  for  bulk  nickel,  as  shown  in  Fig.  22.  It 
melts  at  1950  K,  which  is  higher  than  the  thermodynamic  melting  point  by  222  K.  The  asymptotic 
heterogeneous  melting  point  of  1625  K  is  6  %  lower  than  the  experimental  value. 


Nascent  Nano  Particles 

The  melting  of  nanoscale  particles  is  examined  over  a  size  range  of  1.5-12  nm.  Figure  23  shows 
the  variation  with  temperature  of  the  potential  energy  and  Lindemann  index  of  a  7  nm  aluminum  particle 
consisting  of  10976  atoms.  Melting  occurs  at  810  K.  For  a  particle,  the  surface  acts  as  a  nucleation  site  for 
phase  transition  and  the  melting  front  propagates  to  the  interior  regions  of  the  core.  Hence,  the  phase 
change  is  manifested  by  a  gradual  increase  in  the 
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Figure  23  Potential  energy  and  Lindemann  index  (8)  as  a  function  of  temperature  showing  the 
melting  of  10976-atom  aluminum  (7  nm)  and  nickel  (6  nm)  particles  in  vacuum 
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Figure  24  Effect  of  particle  size  on  melting  of  aluminum  and  nickel  particles 


properties,  as  opposed  to  the  abrupt  changes  found  in  a  bulk  material.  A  similar  trend  is  observed  for  the 
nickel  particle,  except  that  melting  takes  place  at  a  higher  temperature  of  1540  K.  Figure  24  shows  the 
effect  of  particle  size  on  the  melting  of  aluminum  and  nickel  at  nano-scales.  It  also  shows  comparison 
with  the  results  of  previous  MD  simulations  [75,98],  experiments  [99,100],  and  theoretical  studies 
[101,102].  The  theoretical  curve  in  Ref.  [40]  corresponds  to  the  case  in  which  only  the  Ginzburg-Landau 
(GL)  equation  is  employed.  For  nano-sized  particles,  surface  pre-melting  phenomenon  is  observed.  Since 
the  percentage  of  surface  atoms  increases  with  decreasing  particle  size,  a  larger  particle  melts  at  a  higher 
temperature.  Note  that  the  dependence  of  the  melting  point  on  particle  size  becomes  much  weaker  for 
diameters  greater  than  5  nm.  Eckert  et  al.  [99]  synthesized  nanocrystalline  aluminum  powders  by 
mechanical  attrition  in  argon,  hydrogen,  and  oxygen  atmospheres.  A  similar  reduction  was  observed  in 
the  melting  point  with  decreasing  grain  size.  It  is  rather  surprising  that  the  predicted  values  exhibit 
excellent  agreement  with  experimental  data,  since  the  model  under-predicts  the  bulk  melting  point  by  5 
%.  This  may  be  attributed  to  the  fact  that  the  heating  rates  in  the  simulations  are  orders  of  magnitude 
higher  than  those  employed  in  the  experiments  [103].  Furthermore,  the  theoretical  and  experimental 
melting  points  increase  more  gradually  than  those  in  MD  simulations.  This  is  completely  justified,  since 
the  percentage  of  surface  atoms  in  the  particle  follow  a  similar  trend.  Qi  et  al.  [98]  calculated  the  melting 
points  of  nickel  particles  as  a  function  of  particle  size  using  the  quantum-corrected  Sutton-Chen  potential. 
As  can  be  seen,  the  Cleri-Rosato  potential  function  offers  more  accurate  values  of  the  melting  point  of 
nickel  particles. 


Nickel-Coated  Aluminum  Particles 

The  theoretical  model  is  also  employed  to  analyze  the  thermo-chemical  behavior  of  nickel-coated 
aluminum  particles.  Figure  25  shows  the  initial  structure  of  nickel-coated  aluminum  particle  with  a  core 
diameter  of  12  nm  and  shell  thickness  of  1  nm.  A  spherical  nickel  particle  of  known  dimension  is  first 
generated.  A  spherical  void  is  then  created  in  the  interior  of  the  nickel  particle  to  accommodate  the 
aluminum  core.  The  resulting  particle  is  equilibrated  at  300  K  prior  to  the  heating  simulation.  Table  4 
shows  the  values  of  the  core  diameter  ( dc ),  shell  thickness  (<55),  total  number  of  atoms  (A),  and  number  of 
aluminum  atoms  (NAi).  The  core  diameter  varies  in  the 
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Figure  25  Initial  crystal  structure  of  nickel-coated  aluminum  particle  (dp  =  dc+2Ss ) 


range  of  3-12  nm,  and  three  different  shell  thickness  of  0.5,  1,  2,  and  3  nm  are  considered.  The 
aluminum  atomic  fraction  increases  with  increasing  core  size  and  decreasing  shell  thickness. 

Table  4  Configuration  of  nickel  coated  nano-aluminum  particles 


dc ,  nm 

Ss,  nm 

N 

Na, 

xAl 

3 

0.5 

2909 

874 

0.30 

3 

1.0 

7187 

874 

0.12 

3 

2.0 

17440 

874 

0.05 

3 

3.0 

37065 

874 

0.02 

5 

1.0 

15496 

4081 

0.26 

6 

1.0 

20914 

6380 

0.30 

6 

0.5 

12760 

6380 

0.50 

8 

1.0 

39435 

15504 

0.39 

12 

1.0 

103737 

53752 

0.52 

12 

2.0 

169899 

53752 

0.31 

12 

3.0 

254614 

53752 

0.21 
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Baseline  Simulation 

To  facilitate  detailed  discussion  on  the  thermochemical  behavior  of  nickel-coated  aluminum 
particles,  a  particle  with  a  core  diameter  of  12  nm  and  shell  thickness  of  1  nm  is  first  considered.  Figure 
26  shows  the  variation  of  the  Lindemann  index  and  potential  energy  of  the  core  with  temperature.  Melting 
of  the  core  is  characterized  by  the  change  in  the  parameters,  beginning  at  1000  K.  The  predicted  melting 
point  is  higher  than  that  of  a  nascent  particle,  which  is  approximately  880  K.  It  is,  however,  lower  than 
the  structural  melting  point  of  bulk  aluminum,  1060  K.  The  latter  represents  the  upper  limit  for  the  core 
melting  point,  since  it  is  the 


Figure  26  Lindemann  index,  bAh  and  potential  energy  of  the  aluminum  core  as  a  function  of 
temperature  for  a  14  nm  particle  with  a  1  nm  thick  nickel  shell 


temperature  at  which  the  crystal  undergoes  catastrophic  mechanical  failure.  The  melting  point  elevation 
for  a  coated  nano-particle  has  been  observed  in  previous  experimental  [97,104-107]  and  theoretical 
[108,109]  studies.  The  enhanced  thermal  stability  of  an  encapsulated  core  is  attributed  to  the  epitaxial 
core-shell  interface  and/or  to  the  pressure  build  up.  The  former  is  present  in  the  particles  considered  in  the 
present  study,  as  seen  in  Fig.  25.  In  an  uncoated  particle,  the  surface  atoms,  which  have  higher  energy 
than  the  atoms  in  the  interior  region,  can  vibrate  easily.  For  a  coated  particle,  the  shell  restricts  the 
vibration  of  the  interfacial  core  atoms,  which  have  lower  energy  than  those  in  the  interior  region.  As  a 
result,  a  significant  increase  in  the  melting  point  is  observed. 

Figure  27  shows  the  snapshots  of  the  particle  at  300  and  1 100  K.  Melting  causes  a 
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Figure  27  Snapshot  of  the  particle  before  (300  K)  and  after  (1100  K)  core  melting 

significant  change  in  the  crystal  structure  of  the  core.  At  300  K,  a  well-defined  structural  order  is  present, 
but  no  such  order  is  observed  at  1100  K.  The  crystal  structure  of  the  shell  is  negligibly  affected.  Figure 
28  shows  the  variation  of  the  core  radius  with  temperature.  The  sudden  increase 


T  emp  erature.  K 

Figure  28  Core  radius  as  a  function  of  temperature  for  14  nm  particle  with  a  1  nm  thick  Ni  shell 


in  the  core  radius  at  -1000  K  may  be  attributed  to  the  volume  dilation  of  the  core  associated  with  melting. 
Note  that  the  aluminum  core  dissolves  nickel  atoms  at  the  interface.  The  aluminum  atoms  continue  to 
diffuse  into  the  shell,  which  causes  the  core  radius  to  increase  monotonically.  Figure  29  shows  the 
diffusion  processes  occurring  in  the  particle  over  the  temperature  range  of  interest  (300-2400  K). 
Aluminum  atoms  diffuse  into  the  shell  and  nickel  atoms  diffuse  into  the  core,  resulting  in  the  formation  of 
a  homogeneous  alloyed  particle.  The  diffusion  processes  prevail  upon  melting  due  to  the  higher  mobility 
of  the  melted  atoms.  Note  that  complete 
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Figure  29  Snapshots  of  14  nm  particle  showing  melting  and  diffusion  processes 

homogenization  is  not  observed  at  2400  K,  suggesting  that  higher  temperatures  are  required  for  complete 
mixing  of  the  core  and  shell  atoms.  Figure  30  shows  the  variation  of  the  average  potential  energy  of  the 
particle  with  temperature.  The  potential  energy  rises,  attains  a  plateau,  and  then  decreases.  This  trend  is  a 
characteristic  of  the  core-shell  particle  structure,  and  is  not  observed  for  nascent  particles.  The  result 
suggests  the  formation  of  low-energy  species  from  the  following  exothermic  inter-metallic  reaction: 


T  emperature,  K 

Figure  30  Average  potential  energy  of  the  particle  as  a  function  of  temperature  for  core  diameter  of 

12  nm  and  shell  thickness  of  1  nm 
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Ni  +  Al->NiAl. 


(22) 


The  heat  of  formation  of  NiAl  is  62  kJ/mol  at  room  temperature  [110].  The  initial  rise  in  the  potential 
energy  is  caused  by  the  transfer  of  energy  from  the  heat  reservoir  to  the  particle.  The  plateau  represents 
the  stage  at  which  this  energy  supply  is  counterbalanced  by  the  formation  of  the  low-energy  species, 
NiAl.  The  subsequent  decrease  indicates  the  preponderance  of  the  inter-metallic  reactions. 

Note  that  the  interfacial  aluminum  atoms  mix  with  nickel  atoms  even  at  room  temperature.  If  the 
particle  is  heated  externally  to  a  temperature  at  which  the  reactions  begin  to  occur  rapidly,  it  can  ignite, 
provided  that  the  rate  of  chemical  heat  generation  is  greater  than  that  of  heat  loss  to  the  ambient 
environment.  An  energy  balance  is  performed  to  estimate  the  equilibrium  temperature  of  the  particle  upon 
completion  of  the  inter-metallic  reactions  at  adiabatic  conditions.  The  result  is  compared  with  the  value 
obtained  from  the  MD  simulation.  The  thermodynamic  energy  balance  takes  the  form 

H^(T,)  =  HprM(Tj,  (23) 

where  Hreac  is  the  enthalpy  of  the  reactants  calculated  at  an  initial  temperature,  Th  and  Hprod  the  enthalpy 
of  the  products  evaluated  at  the  adiabatic  reaction  temperature,  Tad.  The  initial  temperature  is  taken  as 
1250  K,  which  corresponds  to  the  point  at  which  the  potential  energy  begins  to  decrease.  The  enthalpy  of 
the  reactants,  Hreac ,  is  given  by 
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(24) 


The  enthalpies  of  Al  and  Ni  are  calculated  as  22.85  and  24.82  kJ/mol,  respectively,  and  that  of  the 
melting  of  aluminum  is  taken  as  10.71  kJ/mol.  The  enthalpy  of  the  products  is  expressed  as 
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where  A  is  the  interfacial  area,  V  the  core  volume  and  tt  the  thickness  of  the  interfacial  zone.  The 
interfacial  core  atoms  participate  in  alloying  reactions  prior  to  melting  of  the  core,  thereby  resulting  in  a 
decrease  in  the  energy  content  of  the  particle.  The  fraction  of  the  core  volume  that  has  already  reacted 
depends  on  the  interfacial  area-to-core  volume  ratio  and  the  thickness  of  the  interfacial  zone.  The  latter  is 
approximated  to  be  0.7A  [80].  The  specific  heat  and  enthalpy  of  melting  of  NiAl  are  taken  as  55  J/mol-K 
and  63  kJ/mol,  respectively  [111].  The  equilibrium  temperature  calculated  using  Eqs.  (13)  and  (14)  is 
2388  K.  In  other  words,  the  particle  is  self-heated  from  1250  to  2388  K  due  to  the  heat  release  from  inter- 
metallic  reactions. 

Figure  31  shows  the  variation  of  the  temperature  of  the  particle  with  time  under  adiabatic 
conditions  calculated  using  MD  simulations.  The  initial  position  and  velocities  of  atoms  are  those 
obtained  from  the  heating  simulation  in  an  NPT  ensemble  at  1250  K.  The  particle  is  heated  from  1250  to 
-2300  K  due  to  exothermic  intermetallic  reactions  of  aluminum  and  nickel  atoms  over  a  period  of  100  ps. 
One  of  the  reasons  for  the  discrepancy  is  that  equilibrium  is  not  attained  within  a  time  period  of  -100  ps. 
After  a  time  period  of  50  ps,  the  temperature  increases  very  slowly.  The  particle  snapshots  obtained  in 
NPT  simulations  also  support  this  argument.  It  is  also  important  to  note  that  the  atomic  species  are  no 
longer  expected  to  interact  with  a  ground  electronic  configuration,  but  with  an  excited  one  at  higher 
temperatures.  Therefore,  classical  many-body  potentials  reproduce  interaction  forces  only  with  modest 
accuracy,  and  ab  initio  methods  should  be  applied.  Note  that  the  particle  temperature  increases  by  -1000 
K  over  a  time  period  of -100  ps.  The  resulting  intrinsic  heating  rate  is  -  1013  K/s,  which  is  equivalent  to 
particle  ignition. 
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Figure  31  Temporal  evolution  of  the  particle  temperature  for  core  diameter  of  12  nm  and  1  nm 

thick  shell 


Effect  of  Core  Size 

The  effect  of  core  size  is  studied  in  the  diameter  range  of  3-12  nm  with  a  fixed  shell  thickness  of 
1  nm.  Figure  32  shows  the  variation  of  the  melting  point  of  the  core  with  the  core  diameter.  It  increases 
from  775  to  1000  K  when  the  core  diameter  increases  from  3  to  12  nm.  Such  a  trend  has  also  been 
observed  in  a  previous  study  [78].  In  all  the  cases,  the  presence  of  the  nickel  shell  increases  the  melting 
point  of  the  core.  Figure  33  shows  the  snapshots  of  8  nm  particle  colored  by  thermal  displacement  at 
different  temperatures.  As  can  be  seen,  the  nucleation  of  the  liquid  phase  begins  at  the  core-shell  interface 
and  the  melting  front  propagates  to  the 


Core  Diameter,  nm 

Figure  32  Effect  of  core  diameter  on  the  melting  point  of  the  core  for  a  shell  thickness  of  1  nm 

interior  region  of  the  core.  The  fact  that  heterogeneous  nucleation  at  the  interface  dominates  melting  of 
superheated  nano-particles  has  been  observed  previously  in  MD  simulations  [108].  Figure  34  shows  the 
variation  of  the  core  radius  with  temperature  for  different  core  sizes.  For  a  3  nm  core,  the  diffusion  rate  of 
aluminum  atoms  increases  suddenly  at  1400  K,  a  phenomenon  that  is  not  related  to  the  melting  of  the 
core.  The  nickel  shell  melts  at  1350  K, 
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Figure  33  Particle  snapshots  showing  thermal  displacement  of  core  atoms  at  different  temperatures 

for  8  nm  particle  with  1  nm  shell 
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Figure  34  Normalized  core  radius  as  a  function  temperature  for  shell  thickness  of  1  nm 


which  is  lower  than  the  bulk  melting  point  of  nickel,  1728  K.  Unlike  a  particle,  a  shell  can  begin  to  melt 
both  at  the  surface  and  at  the  core-  shell  interface.  The  nickel  shell  with  a  molten  aluminum  cluster  is  less 
stable  than  a  pure  nickel  particle  of  the  same  size  [97].  For  larger  core  diameters,  the  aluminum  atoms 
begin  to  diffuse  outwards  at  lower  temperatures.  Figure  35  shows  the  variation  of  the  average  potential 
energy  of  the  particle  with  temperature.  The  magnitude  by  which  the  potential  energy  decreases  is 
proportional  to  the  core  size.  In  other  words,  the  smaller  the  core,  the  lower  the  chemical  heat  release  due 
to  the  inter-metallic  reactions.  The  temperature  at  which  the  potential  energy  begins  to  decrease  may  be 
identified  as  the  ignition  temperature.  The  calculated  values  are  in  the  range  of  1280-1480  K,  which  are 
near  the  melting  point  of  the  shell. 


Figure  35  Effect  of  core  diameter  on  the  potential  energy  for  a  shell  thickness  of  1  nm 

Simulations  in  NVE  ensembles  are  performed  to  calculate  the  equilibrium  temperature  of  the 
particle  upon  completion  of  the  inter-metallic  reactions  for  different  core  sizes.  Figure  36  shows  the  effect 
of  core  diameter  on  the  adiabatic  reaction  temperature  of  the  particle.  At  a  fixed  shell  thickness,  the 
reaction  temperature  decreases  with  decreasing  core  size.  For  smaller  cores,  the  temperature  rise  due  to 
the  inter-metallic  reactions  is  relatively  low.  To  understand  the 
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Figure  36  Effect  of  core  diameter  on  adiabatic  reaction  temperature  for  1  nm  thick  shell 

relevant  physiochemical  phenomena,  a  companion  thermodynamic  analysis  is  also  performed  and  the 
result  is  included  in  the  figure.  Figure  37  shows  the  equilibrium  phase  diagram  for  the  Ni-Al  system.  The 
particle  composition  determines  the  inter-metallic  compounds  that  are  formed.  For  3  nm  core,  the 
aluminum  atomic  fraction  is  0.12  (nickel-rich).  The  resulting  product  is  a  mixture  of  Ni3Al  and  Ni.  For 
this  case,  the  chemical  reaction  can  be  expressed  as 

7.2Ni  +  A1  ->  Ni3Al  +  4.2Ni.  (26) 


Figure  37  Ni-Al  phase  diagram  (  s  :  MA13,  d:  Ni2Al3,  J3:  NiAl,  0:  Ni5Al3,  and  a:  Ni3Al)  [112]. 

In  the  present  analysis,  the  specific  heat,  enthalpy  of  fusion  and  heat  of  formation  of  Ni3Al  are 
taken  to  be  130  J/mol-K,  50  kJ/mol,  and  -164  kJ/mol,  respectively  [52].  The  results  support  the  fact  that 
the  adiabatic  reaction  temperature  decreases  with  decreasing  core  size.  This  can  be  attributed  to  the 
following  reasons.  First,  residual  Ni  atoms  remain  in  the  products,  which  do  not  participate  in  the 
reactions.  Second,  although  the  interfacial  zone  thickness  is  held  constant,  the  fraction  of  interfacial  core 
atoms  increases  with  decreasing  particle  size.  As  a  result,  a  greater  percentage  of  core  atoms  have  already 
reacted  for  a  smaller  core.  Such  an  effect  was  also  observed  for  aluminum-coated  nickel  particles  [80]. 
For  aluminum-rich  particle,  inter-metallic  compounds  such  as  NiAl3  are  formed.  The  situation  is 
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typically  encountered  for  aluminum-coated  nickel  particles.  The  ensuing  thermal  behavior  is  expected  to 
be  different. 

Effect  of  Shell  Thickness 

The  effect  of  shell  thickness  is  examined  in  the  range  of  0. 5-3.0  nm.  Figure  18  shows  the  melting 
point  of  the  aluminum  core  as  a  function  of  shell  thickness.  The  core  melting  point  is  not  significantly 
affected  by  variations  in  the  shell  thickness  in  the  size  range  of  concern. 


Figure  38  Effect  of  shell  thickness  on  the  melting  temperature  of  the  12  nm  core 

Calculations  indicate  that  the  core  is  under  compression  (positive  pressure  on  the  order  of  1000  MPa)  and 
the  shell  under  tension.  A  similar  trend  was  also  obtained  for  homogeneously  heated  oxide-coated 
aluminum  particles  [114].  Levitas  [115]  proposed  a  melt  dispersion  mechanism,  which  is  valid  at  very 
high  heating  rates  (>  106  K/s).  According  to  this  theory,  the  shell  is  ruptured  by  the  tensile  stress  exerted 
by  the  core  and  the  molten  aluminum  clusters  are  dispersed  in  the  surrounding  environment.  In  the 
present  study,  melt  dispersion  is  not  observed,  thereby  suggesting  that  the  stresses  are  relaxed.  One 
possible  mechanism  for  such  a  phenomenon  is  the  diffusion  of  aluminum  atoms.  Previous  studies 
[116,117]  indicate  that  the  independence  of  melting  temperature  with  shell  thickness  result  from  stress 
relaxation.  It  is,  thus,  not  surprising  that  the 
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Figure  39  Variation  of  the  core  radius  with  temperature  for  a  core  diameter  of  3  nm  and  shell 

thickness  in  the  range  of  0.5-3  nm 
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melt  dispersion  mechanism  is  not  observed  for  nickel-coated  aluminum  particles.  Figure  39  shows  the 
core  radius  as  a  function  of  temperature  for  a  core  diameter  of  3  nm  and  shell  thickness  of  0.5,  1,  2,  and 
3nm.  Diffusion  is  facilitated  for  thinner  shells.  For  a  0.5  nm  shell,  the  core  radius  increases  significantly 
upon  melting  of  the  core  at  800  K.  The  core  atoms,  however,  begin  to  diffuse  out  only  at  1400,  1500,  and 
1600  K  for  shell  thickness  of  1,  2  and  3  nm,  respectively.  For  a  core  diameter  of  3  nm,  the  results  suggest 
that  the  ignition  temperature  increases  from  800  to  1600  K,  when  the  shell  thickness  increases  from  0.5  to 
3.0  nm.  Figure  40 


Figure  40  Effect  of  shell  thickness  on  the  melting  temperature  of  the  nickel  shell  for  core  diameters 

of  3  nm. 

shows  the  effect  of  shell  thickness  on  the  melting  point  of  the  shell  for  core  diameter  of  3  nm.  The  shell 
melting  point  increases  from  a  value  as  low  as  1320  to  1580  K,  when  the  thickness  increases  from  1  to  3 
nm.  In  all  cases,  the  melting  point  of  the  shell  is  lower  than  that  corresponding  bulk  value  of  1728  K.  It  is 
obvious  that  shell  melting  plays  an  important  role  in  dictating  the  diffusion  and  ignition  in  nickel-coated 
aluminum  particles. 

The  case  of  a  particle  covered  by  a  very  thin  shell  is  of  interest,  since  a  previous  study  indicates 
the  possibility  of  ignition  near  room  temperature  [81].  Furthermore,  it  is  possible  to  create  a  particle  with 
a  reasonably  large  core  size  and  a  thin  shell  that  features  an  equal  number 
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Figure  41  Diffusion  coefficient  of  aluminum  in  nickel  for  a  7  nm  particle  with  a  0.5  nm  shell 
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of  A1  and  Ni  atoms.  As  a  consequence,  the  energetic  of  inter-metallic  reactions  are  easily  tractable.  In  the 
present  study,  a  7  nm  particle  with  a  0.5  nm  shell  is  considered.  The  information  about  the  diffusion 
coefficient  of  aluminum  atoms  in  a  nickel  shell  is  useful,  especially  from  the  standpoint  of  development 
of  macro-scale  ignition  and  combustion  models.  Figure  41  shows  the  effect  of  temperature  on  the 
diffusion  coefficient  of  aluminum  atoms  in  a  nickel  shell.  The  diffusionprocess  is  very  slow  at 
temperatures  lower  than  the  melting  point  of  the  core.  A  sudden  increase  in  the  diffusion  coefficient  takes 
place  after  melting  of  the  core.  It  exhibits  a  temperature  dependence  of  the  form 


f 

D  =  D0  exp 

v 


(27) 


where  EA  is  the  activation  energy  and  D0  is  the  pre-exponential  factor.  The  curve-fit  indicates  a  value  of 
34.7  kJ/mol  for  the  activation  energy  and  8.19xl0"8  m2/s  for  the  pre-exponential  factor.  The  obtained 
diffusion  coefficients  can  be  used  to  study  combustion  of  nickel-coated  aluminum  particles.  Figure  42 
shows  the  temporal  evolution  of  the  particle  in  adiabatic  condition.  The  temperature  increases  suddenly 
upon  melting  of  the  core.  The  inter-metallic  reactions  result  in  particle  heating  from  300  to  1420  K.  The 
results  substantiate  the  fact  that  an  aluminum  particle  coated  by  a  thin  nickel  shell  can  be  pyrophoric  even 
in  an  inert  environment. 


Time,  ns 

Figure  42  Temperature  as  a  function  of  time  for  a  7  nm  particle  with  a  0.5  nm  thick  shell  under 

adiabatic  condition 

The  thermo-chemical  behavior  of  nickel-coated  aluminum  particles  in  the  size  range  of  4-18  nm 
was  studied  using  molecular  dynamics  simulations.  The  analysis  was  carried  out  in  isothermal-isobaric 
and  isochoric-isoenergetic  ensembles  using  an  embedded  atom  method.  Emphasis  was  placed  on 
analyzing  the  melting  points  of  the  core  and  shell,  diffusion  of  atoms,  and  inter-metallic  reactions.  The 
aluminum  core  melted  at  a  temperature  greater  than  the  melting  point  of  a  nascent  aluminum  particle  due 
to  the  cage-like  mechanical  constraint  imposed  by  the  nickel  shell.  The  melting  point  of  the  aluminum 
core  increased  from  775  K  at  3  nm  to  1000  K  at  12  nm.  It  was  not  significantly  affected  by  the  variations 
in  the  shell  thickness  in  the  range  of  1-3  nm.  The  melting  point  of  the  shell  increased  with  increasing 
thickness,  from  a  value  as  low  as  1 100  K  at  1  nm  to  1580  K  at  3  nm.  Melting  was  followed  by  diffusion 
of  atoms  and  energy  release  due  to  inter-metallic  reactions,  which  resulted  in  self-heating  of  the  particle. 
The  diffusion  coefficient  of  aluminum  atoms  in  the  nickel  shell  had  an  exponential  dependence  on 
temperature,  with  an  activation  energy  of  34.7  kJ/mol.  The  adiabatic  reaction  temperature  of  the  particle 
increased  from  1650  to  2338  K,  when  the  core  diameter  increased  from  3  to  8  nm.  The  calculated  values 
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exhibited  reasonably  good  agreement  with  predictions  from  a  thermodynamic  energy  balance  analysis. 
The  results  demonstrated  that  nickel-coated  aluminum  particles  can  ignite  even  in  inert  environments  (or 
in  the  absence  of  an  oxidizing  gas).  Future  work  will  focus  on  the  energetics  of  aluminum-coated  nickel 
particles  and  the  comparison  with  the  results  of  nickel-coated  aluminum  particles  obtained  from  the 
present  study.  Furthemore,  a  multi-scale  model  will  be  developed  to  study  the  flame  propagation  of 
nickel-coated  aluminum  particles.  Results  from  MD  simulations  will  be  used  to  characterize  the  burning 
behavior  of  individual  particles. 

B.  Combustion  of  Nano- Aluminum  Particles  and  Liquid  Water  Mixtures  (Georgia 

Tech  with  Penn  State) 

The  combustion  of  aluminum  particles  in  water  is  of  relevance  to  many  propulsion  and  energy- 
conversion  applications.  In  metalized  composite  solid  propellants,  aluminum  particles  typically  react  with 
the  combustion  products  of  the  polymeric  binder  and  ammonium  perchlorate,  of  which  water  vapor  (H20) 
and  carbon  dioxide  (C02)  are  two  major  species  [118].  The  problem  is  of  particular  interest  for 
underwater  propulsion,  since  the  oxidizer  (water)  could  be  supplied  from  the  environment  [119].  The 
combustion  of  aluminum  particles  in  water  has  also  been  studied  in  the  contexts  of  hydrogen  generation 
[120],  nuclear  reactor  and  industrial  manufacturing  safety  [121],  and  underwater  explosives  [122].  The 
frozen  mixture  of  aluminum  particles  and  water  (ALICE)  is  also  under  consideration  for  several  energetic 
applications,  due  to  its  structural  integrity  ([123],  [124]). 


Figure  43  Effect  of  pressure  on  aluminum  vaporization  temperature  and  adiabatic  flame 

temperatures  for  various  oxidizers 

The  ignition  and  combustion  characteristics  of  aluminum  particles  in  water  are  different  from 
those  in  oxygenated  environments.  Figure  43  shows  the  effect  of  pressure  on  the  aluminum  vaporization 
temperature  ( Tv>Ai )  and  adiabatic  flame  temperature  (7})  of  aluminum- water  mixtures  with  different  particle 
sizes  and  water  in  various  thermodynamic  states.  The  particles  are  assumed  to  be  passivated,  with  an 
oxide  shell  thickness  of  3  nm.  The  flame  temperature  is  lower  than  the  vaporization  temperature  for 
pressures  over  a  “cut-off’  value.  For  micron-sized  particles,  this  occurs  in  the  range  of  -1-4  atm, 
depending  on  the  thermodynamic  state  of  water.  For  nano-sized  particles,  the  cut-off  pressure  is  as  low  as 
-0.2  atm.  This  can  be  attributed  to  the  fact  that  the  inert  oxide  layer  constitutes  a  greater  portion  of  the 
particle  mass  at  nano-scales;  a  50  nm  aluminum  particle,  for  example,  contains  32  wt.  %  oxide  [122].  It 
should  be  pointed  out  that  the  actual  flame  temperature  can  be  lower  than  the  theoretical  value,  due  to  the 
effects  of  heat  losses  and  incomplete  combustion.  As  a  result,  in  most  practical  cases,  nano-aluminum 
does  not  vaporize  and  heterogeneous  chemical  reactions  occur  at  the  particle  surface. 
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Figure  44  shows  the  effect  of  particle  size  on  the  ignition  temperature  of  aluminum  particles  in 
water- vapor  containing  environments.  Oxidation  was  studied  by  various  methods  including  a 
thermogravimetric  analyzer  [125],  electrical  heating  [126],  shock  tube  [127],  hydrogen-oxygen-argon 
burner  [128],  and  arc  burner  [129].  Both  spherical  particles  ([125],  [127-129])  and  wire  samples  [126] 
were  considered.  In  the  thermo-gravimetric  experiments  [125],  the  particles  were  heated  at  relatively  low 
rates  in  the  range  of  1-20  K/min.  The  oxidizing  gas  consisted  of  27  %  H20  and  73  %  Ar.  The  reactions 
were  observed  at  much  lower  temperatures  in  water  vapor  than  in  oxygen.  Specifically,  a  stepwise  weight 
change  was  observed  at  the  melting  point  of  aluminum  (660  °C)  and  the  particle  is  completely  oxidized  at 
-1000  °C.  Note  that  the  experimental  conditions  are  different  in  the  referenced  studies.  For  example,  the 
heating  rates  in  shock  tube  experiments  [127]  are  much  higher  (-106  K/s)  and  the  oxidizer  concentration 
is  different  from  those  employed  in  other  studies.  It  is  also  worthwhile  to  mention  that  the  measured 
ignition  temperature  in  burner  ([128],  [129])  and  shock  tube  [127]  experiments  corresponds  to  the 
temperature  of  the  gas.  The  ignition  temperature  depends  on  various  parameters  such  as  the  particle 
characteristics,  oxidizer  concentration,  and  heating  rate.  For  these  reasons,  the  experimental  data  cannot 
be  quantitatively  compared  directly.  Nevertheless,  it  is  true  that  the  ignition  temperature  increases  with 
increasing 
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Figure  44  Effect  of  particle  size  on  ignition  temperature  of  aluminum  particles  in  water  vapor 


particle  size  in  the  range  of  10  nm-100  pm.  Even  for  larger  particles,  it  remains  lower  than  the  melting 
point  of  the  oxide  layer,  2350  K.  This  suggests  that  the  ignition  temperatures  of  aluminum  particles  in 
water  vapor  are  lower  than  those  in  oxygen.  The  phenomenon  may  be  attributed  to  the  stabilization  of  y- 
oxide  polymorph  [125]  and/or  formation  of  a  weaker  hydroxide  layer  on  the  particle  surface  [126]. 

Many  of  the  desirable  traits  of  nano-sized  particles  are  due  to  the  presence  of  a  large  percentage 
of  atoms  on  the  surface  [130].  The  percentage  of  atoms  on  the  surface  layer  of  an  aluminum  particle 
increases  from  2  to  92  %,  when  the  particle  size  decreases  from  100  to  1  nm.  The  surface  atoms  have 
lower  coordination  numbers  and  higher  energies  than  the  atoms  in  the  interior  region  of  the  particle.  As  a 
result,  the  thermophysical  properties  of  nano-sized  particles  are  significantly  different  from  their 
corresponding  bulk  values.  The  melting  temperature  of  nano-aluminum  particles  decreases  from  933  to 
473  K  when  the  particle  size  decreases  from  10  to  2  nm  ([131],  [1132]).  Similarly,  the  melting  point  of 
the  oxide  layer  can  be  lower  than  the  bulk  value  of  2350  K;  for  a  shell  thickness  of  2.5  nm  the  melting 
point  of  the  oxide  layer  is  1313K[133].  The  ignition  temperatures  and  burning  times  of  nano-aluminum 
particles  are  also  lower  than  their  micron-sized  counterparts  [134].  Significant  enhancement  in  the 
burning  characteristics  is,  thus,  expected,  when  nano-sized  particles  are  used  to  formulate  energetic 
materials  [135]. 


51 


The  combustion  of  nano-sized  aluminum  particles  and  water  has  been  studied  experimentally  for 
a  relatively  wide  range  of  pressure  and  particle  size  ([122],  [136],  [137]).  Results  from  earlier  studies 
([136], [137])  suggested  that  the  presence  of  a  gelling  agent  such  as  polyacrylamide  is  necessary  to 
achieve  self-sustained  deflagration.  Risha  et  al.  ([122], [138]),  however,  demonstrated  self-deflagration  of 
aluminum-water  mixtures  in  a  constant  volume  optical  pressure  vessel.  The  particle  diameters  covered  a 
range  of  38-130  nm  and  the  pressure  range  of  interest  was  0.1-10  MPa.  The  measured  burning  rates  were 
found  to  be  inversely  proportional  to  particle  size  and  exhibited  a  pressure  dependence  of  the  form  rb  = 
apm ,  with  the  exponent  in  the  range  of  0.27-0.47.  Both  transport  and  chemistry  were  speculated  to  affect 
the  burning  behavior,  but  key  mechanisms  and  parameters  are  yet  to  be  explored  [122].  In  the  present 
study,  a  theoretical  model  is  developed  to  study  the  flame  propagation  of  a  quasi-homogeneous  mixture  of 
nano-sized  aluminum  particles  and  liquid  water.  Special  attention  is  placed  on  the  effects  of  particle  size 
and  pressure  on  the  flame  structure  and  burning  property.  Results  suggest  that  the  combustion  of  nano¬ 
aluminum/water  mixtures  is  controlled  by  mass  diffusion  through  the  oxide  layer  of  the  particles. 
Reasonably  good  agreement  with  experimental  data  is  achieved,  demonstrating  the  validity  of  the 
proposed  model. 

Experiment 

The  burning  rates  of  nano-aluminum/water  mixtures  were  obtained  in  an  argon  environment 
using  a  constant  volume  vessel,  as  shown  schematically  in  Fig.  45.  The  chamber,  made  of  stainless  steel, 
is  equipped  with  four  optical  viewing  ports,  each  having  a  15.2  X  2.54  cm  field  of 
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Figure  45  Schematic  of  constant-pressure  strand  burner  with  optical  access 

view.  The  61 -cm  long  chamber  has  an  inner  diameter  of  22  cm  and  a  total  free  volume  of  23  liters;  the 
relatively  large  volume  minimizes  the  pressure  variation  caused  by  the  generation  of  gaseous  combustion 
products.  The  base  plate  has  six  feed-through  ports  to  provide  pathways  for  electrical  signals  and  gas  lines 
into  the  chamber.  Nano-aluminum  particles  were  obtained  from  Technanogy  and  Nanotechnology.  They 
were  mixed  with  stoichiometric  amount  of  distilled  water  in  a  sealed  plastic  bag.  Stoichiometry  was 
calculated  based  on  the  active  aluminum  content  in  the  particle.  No  gelling  agent  was  employed.  The 
mixture  was  then  packed  into  a  quartz  glass  tube  (1  cm  OD,  0.8  cm  ID,  7.5  cm  long).  Ignition  was 
achieved  using  a  small  -1/8  in.  thick  propellant  booster  made  of  a  homogenous  double-base  propellant 
(NOSOL  363)  initiated  by  a  resistance-heated  nichrome  wire  threaded  through  the  booster.  The  temporal 
evolution  of  the  regressing  luminous  front  was  tracked  and  recorded  using  video  equipment.  Figure  46 
shows  captured  images  of  stoichiometric  mixture  containing  80  nm  particles  burning  at  a  pressure  of  5.8 
MPa.  The  flame  front  propagates  through  the  packed  mixture  at  a  constant  velocity.  The  measured 
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position-time  curve  was  used  to  determine  the  burning  rate  of  the  strand.  More  details  of  the  experimental 
set-up  can  be  found  in  [122]. 


Figure  46  Captured  images  of  stoichiometric  aluminum-water  mixture  containing  80  nm  particles 

burning  at  a  pressure  of  5.8  MPa 


Theoretical  Framework 

The  analysis  considers  one-dimensional,  isobaric,  and  planar  flame  propagation  in  a  chemically 
reacting  system  consisting  of  passivated  nano-aluminum  particles  and  liquid  water.  Note  that  the  present 
study  deals  with  the  downward  propagating  flame.  For  other  conditions,  the  released  hydrogen  gas  is 
likely  to  affect  the  heat  transfer  to  the  unburned  mixture.  The  aluminum-water  mixture  is  a  viscous  paste 
with  dense  particle  loading  and  significant  particle-particle  interactions,  as  opposed  to  dust  clouds  that 
feature  dilute  particle  concentrations  ([134],  [139]).  Furthermore,  the  gaseous  oxidizer  (water  vapor)  is 
produced  in  situ  by  vaporization,  instead  of  being  supplied  externally.  The  particles  are  assumed  to  be 
uniformly  sized  and  their  agglomeration  is  neglected.  Flame  propagation,  which  takes  place  through 
thermal  conduction,  is  accompanied  by  a  variety  of  physicochemical  processes,  including  water 
vaporization,  chemical  reactions,  and  mass,  momentum  and  energy  exchanges  between  the  fluid  and 
particle  phases.  The  system  is  approximated  to  be  pseudo-homogenous,  so  that  the  particles  and 
surrounding  fluid  are  in  thermal  equilibrium  locally.  Water  vaporization  occurs  on  an  infinitesimally  thin 
plane.  Figure  47  shows  the  multi-zone  flame  structure  considered  in  the  present  study.  The  entire  spatial 
domain  is  divided  into  three  zones  to  demarcate  the  regions  in  which  phase  transition  and  chemical 
reactions  occur.  The  initial  temperature  of  the  mixture  is  298  K.  Water  undergoes  a  thermodynamic 
phase  transition  at  the  vaporization  front,  x  =  -Sy,  where  the  local  temperature  reaches  the  vaporization 
point,  Tv.  The  particles  start  to  bum  once  the  ignition  temperature,  Tigm  is  attained.  Chemical  reactions  are 
neglected  in  the  preheat  zones. 

The  thermal  conductivity  of  the  mixture,  2m,  depends  on  the  loading  density  and  thermal 
conductivity  of  the  particles.  The  following  correlation  provides  the  best  fit  to  experimental  data  for  a 
wide  range  of  particle  volume  fractions  [140]: 


K  =  K  exP 


^-1.50/ 


(28) 


where  X  stands  for  the  thermal  conductivity  and  0  the  volume  fraction.  The  subscripts  m,p,  and 


Al-H2Om  Al-H20(g) 
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Figure  47  Physical  model  and  multi-zone  flame  structure  (O  Al,  O  AI2O3) 


/ refer  to  mixture,  particle,  and  fluid,  respectively.  Equation  (28)  assumes  that  the  fluid  is  non-conducting 
when  compared  with  the  particle.  This  is  a  reasonable  approximation,  since  the  former  value  is  two  or 
three  orders  of  magnitude  lower  than  the  latter.  Furthermore,  it  is  reassuring  to  note  that  the  thermal 
conductivity  predicted  by  the  correlation  agrees  with  the  experimental  data  for  a  wide  range  of  particle 
volume  fractions.  The  oxide  layer  hampers  the  flame  propagation,  due  to  its  inefficiency  as  a  conductor  of 
heat.  The  overall  thermal  conductivity  of  the  passivated  aluminum  particle,  Xp  is  calculated  as  follows 
[141]: 


A Al R 
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b  =  2(R-r)Aox+2rAM. 


(29) 


(30) 


Here  r  is  the  radius  of  the  aluminum  core  and  R  the  outer  radius  of  the  particle.  The  subscripts  Al  and  ox 
refer  to  aluminum  and  oxide,  respectively.  The  thermophysical  properties  of  water  and  hydrogen  are 
given  in  Ref.  [142],  while  those  of  aluminum  and  its  oxide  are  taken  from  Refs.  [143-145]  and  Ref.  [146], 
respectively.  Table  5  summarizes  the  property  data  at  a  baseline  pressure  of  3.65  MPa.  All  properties  are 
evaluated  at  an  average  temperature  in  each  zone.  The  enthalpy  of  vaporization  and  boiling  point  of  water 
are  calculated  as  a  function  of  pressure.  A  continuum- 
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Table  5  Thermophysical  properties  of  different  species  in  three  zones  at  baseline  pressure  of  3.65 

MPa* 


Species 

Thermal  conductivity, 
W/m-K 

Specific  heat, 
kJ/kg-K 

Density, 

kg/m3 

W 

V 

G 

W 

V 

G 

W 

V 

G 

Aluminum 

239 

95 

143 

0.954 

1.260 

1.176 

2700 

Aluminum  Oxide 

22.20 

9.50 

6.40 

0.930 

1.200 

1.303 

4000 

Liquid  water 

0.70 

- 

4.400 

- 

995 

- 

- 

Water  Vapor 

- 

0.07 

0.12 

- 

2.330 

2.580 

- 

8.58 

5.17 

Hydrogen 

- 

- 

0.80 

- 

- 

16.150 

- 

- 

0.56 

*  Tv  =  519  K,  Tign  =  1360  K,  Tf  =  1800  Katp  =  3.65  MPa 
W:  liquid  water  zone;  V:  water  vapor  zone;  G:  gaseous  reaction  zone 

regime  heat  transfer  model  is  considered,  since  the  present  study  deals  with  pressures  representative  of 
those  in  practical  propulsion  systems.  Free  molecular  effects  are  neglected  without  significant  reduction 
in  accuracy  [123].  Mass  and  energy  balances  are  enforced  for  a  differential  element  in  each  zone.  The 
resulting  conservation  equations  are  solved  to  obtain  the  burning  rate  and  temperature  distribution.  The 
formulation  is  developed  based  on  a  coordinate  system  attached  to  the  propagating  flame. 

Energy  balance  for  the  liquid  water  zone 

The  liquid  water  zone  encompasses  the  region  between  the  far  field,  *  =  -  oo,  and  the  vaporization 
front,  x  =  ~  Sv.  The  energy  equation  takes  the  form 

(PaiCp,ai®AI  +  Pox^p,ox^ox  +  PlSpA  )  rb  ~  ^ m,W  ^ 2  ’  ^  ^ 


subject  to  the  boundary  conditions: 


T  -T  'T  =  T  ■ 

x-»-00  u  5  x=-S.  v’ 


(32) 


where  p  is  the  density,  Cp  the  specific  heat,  rb  the  burning  rate,  T  the  temperature,  x  the  space  coordinate, 
and  Sv  the  thickness  of  the  vapor  zone.  The  subscripts  W,  u ,  v,  ox,  and  Iw  refer  to  the  liquid  water  zone, 
unbumed  state,  vaporization,  oxide,  and  liquid  water,  respectively.  An  analytical  solution  to  Eq.  (31)  is 
obtained  for  the  temperature  profile. 

T  =  Tu+(Tv-Tu)QMkw(X  +  Sv)}’  (33) 

where  £w  is  the  ratio  of  the  burning  rate  to  the  thermal  diffusivity,  defined  as 

~rb  /  {^m,w  /  (pCp  )r  j  ?  (34) 
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where  ( pCp )  denotes  the  volume-averaged  product  of  the  density  and  specific  heat  of  the  mixture.  The 
temperature  profile  depends  on  the  burning  rate  and  thickness  of  the  water  vapor  zone,  both  of  which  are 
not  known  a  priori. 

Energy  and  mass  balance  for  the  water  vapor  zone 

Water  vapor  generated  at  the  vaporization  front,  x  =  -^  flows  through  the  interstitial  space  between 
particles.  The  conservation  of  mass  of  water  is  enforced  to  determine  the  velocity  of  water  vapor  at  the 
vaporization  front 

Plvfb  Pwv^wv*  (35) 

Here,  v  is  the  gas  velocity.  The  subscript  wv  denotes  water  vapor.  Equation  (35)  can  be  used  to  express 
the  energy  equation  in  the  following  form: 


Al^ pM^Al  +  PoxCp^x  +  P,WCp^lW)rb  ^  ~ 


w  7  2 


dx 


(36) 


subject  to  the  interfacial  conditions: 


o  -  dT 

x  =  -5v-K-r 

ax 

*  =  0:T  =  Tign, 


_  .  dT 

m  i 

ax 


+  hfe(E>iwPh/b  ’ 


w 


(37) 


where  hfg  is  the  enthalpy  of  water  vaporization  and  Tign  the  ignition  temperature  of  nano-aluminum 
particles.  The  subscript  V  refers  to  the  water  vapor  zone.  The  thickness  of  this  zone  is  obtained  by 
performing  the  heat-flux  balance  at  x  =  -Sv: 

U  Krkr(T,„-T,)  M 


'),  =  T—  l0tr'  ■ 

kv 


1  + 


IV 


^m,W^W  (T  Tu  )  +  hfg Piw®  iwrb 


(38) 


The  temperatures  at  the  interfacial  boundaries  are  matched  to  provide  a  closed- form  solution  to  the  energy 
equation: 


T 


(39) 


Energy  balance  for  the  reaction  zone 

The  stoichiometric  reaction  of  aluminum  particles  with  water  vapor  is  given  by 

2AI  +  3H20  — >  Al203  +  3 H2 .  (40) 


The  properties  are  calculated  by  averaging  their  respective  values  of  the  reactant  and  product  species.  The 
energy  equation  can  be  expressed  as  [134] 


f 

V 


m,G 


d  T  ,  p„Qr 

r  2  +  9 

dx  Tb 


(41) 


where  Qr  is  the  chemical  energy  release  per  unit  mass  of  the  mixture,  and  ib  the  particle  burning 

time.  The  subscripts  G,  m  and  i  refer  to  the  reaction  zone,  mixture,  and  species  /,  respectively.  To 

facilitate  the  analysis,  the  temperature  and  spatial  coordinate  are  normalized  as  follows: 
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X 


(42) 


where  z0  is  the  reference  time  scale  defined  as  the  particle  burning  time  at  a  reference  temperature.  The 
location  y  =  0  is  the  ignition  point  of  particles.  Substituting  the  normalized  variables  defined  in  Eq.  (42) 
into  Eq.  (41),  the  non-dimensional  form  of  the  energy  equation  is  obtained 

d  0  ?  dO  ?  /  ^  Tn 


dy  dy 


Here  k  =  rb  *  y  r0  /  amG  is  the  normalized  burning  rate,  with  a  being  the  thermal  diffusivity.  The 
normalized  heat-release,  //,  is  written  as 

-a  rp  O 

m.Gr'  m^r 


^m,G  {^ign  ) 


The  energy  equation  gives  only  a  partial  description  of  the  underlying  physicochemical  phenomena  in  the 
reaction  zone.  In  particular,  an  equation  for  the  consumption  of  the  particle  mass  is  also  needed. 

Particle  mass  consumption  in  the  reaction  zone 

The  composition  of  the  particle  changes  during  the  course  of  its  reaction  with  water  vapor.  The 
aluminum  content  decreases  progressively  and  a  spherical  oxide  particle  forms  after  complete  oxidation. 
To  characterize  the  combustion  of  individual  particles,  an  equation  for  the  particle  mass  consumption  is 
considered 


where  Mp  is  the  particle  mass.  The  subscript  0  refers  to  the  initial  state.  Equation  (45)  can  be  re-written  in 
terms  of  the  particle  size 


where  dp  is  the  particle  diameter.  For  consistency,  the  particle  diameter  is  normalized  as  follows: 


Substituting  the  normalized  variable  defined  in  Eq.  (47)  into  Eq.  (46),  the  non-dimensional  form  of  the 
particle  mass  consumption  equation  is  obtained 
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The  following  boundary  conditions  are  specified  to  close  the  formulation: 


d0_ 

dy 


,  d£ 

m  / 

g  dy 


> ii=1 


y  =  l:^-  =  0,  rj  =  0. 
dy 


(49) 


The  mass  and  energy  balance  equations  in  the  reaction  zone  are  solved  numerically,  with  the  burning  rate 
treated  as  the  eigenvalue.  A  shooting  technique  is  employed  to  find  the  solution;  the  Newton-Raphson 
iteration  method  is  used  [147].  Numerical  integration  is  achieved  by  means  of  the  Rosenbrock  method 
[147]. 


Heat  release 

The  actual  heat  release  from  particle  burning  is  lower  than  its  theoretical  counterpart  due  to  incomplete 
combustion.  The  combustion  efficiency  of  nano-aluminum/water  mixtures  is  in  the  range  of  80-100%, 
depending  on  the  pressure  and  particle  size  [138].  Heat  loss  to  the  environment  also  occurs  through 
thermal  conduction  and  radiation.  To  incorporate  these  effects  into  the  model,  the  normalized  heat  release 
is  calculated  based  on  the  actual  flame  temperature.  A  simplified  expression  for  the  normalized  heat 
release  is  obtained  by  integrating  the  energy  equation,  Eq.  (16),  and  imposing  the  boundary  conditions 
specified  in  Eq.  (49) 


where  the  subscript  / denotes  the  flame.  Diakov  et  al.  [148]  studied  the  flame  propagation  of  aluminum- 
water  mixtures  in  a  stainless  steel  chamber  equipped  with  thermocouples  at  a  pressure  of  1  atm.  The 
particle  size  is  100  nm  and  the  oxide  layer  thickness  is  1.84  nm.  The  measured  flame  temperature  is  1800 
K,  which  is  significantly  lower  than  the  theoretical  value  of  2790  K  for  a  particle  size  of  100  nm.  The 
combustion  efficiency  is  measured  to  be  87%.  A  similar  scenario  was  observed  for  magnesium-water 
mixtures  [149].  The  combustion  efficiency  increases  with  decreasing  particle  size  [138],  although  the 
adiabatic  flame  temperature  is  lower  for  smaller  particles.  For  simplicity,  a  flame  temperature  of  1800  K 
is  used  to  calculate  the  heat  release  for  all  particle  sizes.  Results  of  the  sensitivity  analysis  indicate  that 
the  burning  rate  increases  modestly  with  increasing  flame  temperature.  At  a  pressure  of  3.65  MPa  and 
particle  size  of  38  nm,  the  burning  rate  increases  from  5.59  to  6.85  cm/s  when  the  flame  temperature 
increases  from  1800  to  2300  K. 

Ignition  temperature  and  burning  time  of  particles 

The  present  model  requires,  as  input  parameters,  the  ignition  temperature  and  burning  time  of 
nano-aluminum  particles.  The  ignition  temperatures  of  Technanogy  aluminum  particles  (24-192  nm)  in 
water  vapor  vary  between  1325  and  1360  K  [128],  as  shown  in  Fig.  44.  For  convenience,  we  assume  the 
ignition  temperature  to  be  1360  K. 

The  combustion  of  aluminum  particles  involves  mass  diffusion  through  the  gas-phase  mixture 
and  oxide  layer,  and  chemical  reactions  between  the  aluminum  atoms  and  oxidizer  molecules.  For 
micron-sized  and  larger  particles,  combustion  is  controlled  by  species  diffusion  through  the  gas-phase 
product  mixture;  the  burning  time  follows  dpLS- law  and  exhibits  weak  dependencies  on  temperature  and 
pressure  [134].  For  particles  smaller  than  1  pm,  the  pressure  and  temperature  of  the  ambient  gas 
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significantly  influence  the  burning  time  [150],  whereas  the  particle  size  exerts  only  a  weak  effect  [134]. 
The  latter  trend  may  be  attributed  to  the  transition  from  the  continuum  to  the  free-molecular  heat  transfer 
regime  [34]  and/or  sintering  and  agglomeration  of  particles.  The  present  study  deals  with  the  conditions  in 
the  continuum  regime,  in  which  the  burning  rate  is  more  size  dependent.  For  nano-sized  particles,  surface 
tension  leads  to  strong  adherence  of  the  oxide  layer  to  the  particle  [152].  This  phenomenon  was  also 
observed  in  molecular  dynamics  simulations  [133]. 

The  diffusion  resistance  provided  by  the  oxide  layer  is  several  orders  of  magnitude  greater  than 
that  of  the  gas-phase  mixture.  As  a  result,  species  diffusion  through  the  oxide  layer  is  the  rate-controlling 
process  [153].  Park  et  al.  [154]  studied  the  oxidation  of  nano-aluminum  particles  in  air  using  single 
particle  mass  spectrometry  for  temperatures  up  to  1373  K.  The  measured  oxidation  rates  exhibited 
reasonably  good  agreement  with  the  results  of  the  diffusion-controlled  combustion  model.  The 
kinetically-controlled  model,  however,  failed  to  capture  the  observed  trends.  These  results  support  the 
hypothesis  that  the  combustion  of  nano-sized  aluminum  particles  is  controlled  by  species  diffusion 
processes  rather  than  by  chemical  kinetics.  To  estimate  the  particle  burning  time,  therefore,  the  diffusion 
coefficients  of  the  reacting  species  in  the  oxide  layer  are  required.  They  are,  however,  poorly  known,  and 
only  with  uncertainties  of  several  orders  of  magnitude  [155].  In  the  present  study,  the  reaction  time  scale 
is  approximated  to  be  the  particle  burning  time.  This  is  physically  justified  because  the  energy  release 
rate  is  dictated  by  the  reactivity  of  the  individual  particles  [134].  Theory  suggests  that  the  characteristic 
time  scale  for  species  diffusion  through  the  oxide  layer  is  inversely  proportional  to  pressure  and  bears  a 
^/-relationship  ([153],  [154]).  The  observed  particle  size  effect  on  the  burning  rate  is  obtained  only  if 
such  a  relationship  is  employed.  A  dp-  law  is,  thus,  adopted  for  the  burning  time.  The  reference  time  scale 
is  taken  to  be  the  burning  time  of  a  24  nm  aluminum  particle  and  the  pressure  exponent  in  the  burning 
time  relationship  varies  between  -0.3  and  -1.0,  depending  on  the  ambient  temperature  ([128],  [148]): 


*b  = 


c\jax  exp (bxT )  +  a2  exp ( b2T )] d2p 


(51) 


m  =  a3  exp (&3r)  +  a4  exp (b4T^j, 


(52) 


where  dp  is  the  particle  diameter  in  nm  and  p  the  ambient  pressure  in  atm.  The  constants  are  given  in 
Table  6.  The  theoretical  studies  on  oxidation  of  nano-aluminum  particles  also  indicate  that  the  diameter 
exponent  in  the  burning  time  relationship  is  in  the  range  1. 6-2.0  [155]. 


Table  6  Constants  in  burning  time  expression 


Constant 

Value 

c 

1.736xl0'3 

ai 

204.650 

bi 

-9.848xl0"3 

a2 

1.842xl0'4 

b2 

3.461xl0'5 

a3 

7.075 

b3 

-1.905xl0'3 

a4 

4.023xl0_1 

b4 

-3. 120x1  O'4 
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Analytical  Model  of  the  Mixture  Burning  Rate 

With  judicious  simplifications,  it  is  possible  to  obtain  an  analytical  expression  for  the  burning  rate 
[149].  The  specific  heat  capacities  of  different  species  and  thermal  conductivities  of  the  mixture  in  the 
three  zones  (see  Fig.  47)  are  assumed  to  be  equal.  The  energy  equation  can  be  written  as 

£-*£.  (53) 

where  &  is  the  volume  fraction  of  water.  For  an  aluminum- water  mixture,  &  takes  the  following  form: 

Pn 

0  =  — ^ — .  (54) 

Pp  +  Plw 

Equations  (38)  and  (39)  are  combined  to  provide  an  expression  for  the  heat  flux  at  x  =  0  in  the  preheat 
zone 
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,  PlwPp 
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-C  (t  )+  PlwPp  h 

p\  w  «)  n  „  - 


Plw  +  Pp 


fg 


(55) 


The  heat  flux  at  x=0  in  the  reaction  zone  is  calculated  by  taking  the  spatial  derivative  of  the  analytical 
solution  to  Eq.  (41) 


xdT_  XQr 
dx  G  TbrbCp  ’ 


(56) 


By  matching  the  two  heat  fluxes  at  x  =  0,  an  analytical  expression  for  the  burning  rate  is  obtained 


r  -  LA _ M _ 1  (57) 

‘  V-c*  2CPK-T.)+K  p 

where  pm  is  the  density  of  the  unbumed  mixture.  The  parameters  that  dictate  flame  propagation  are  the 
thermal  diffusivity  of  the  mixture,  enthalpy  of  reaction,  ignition  temperature,  and  burning  time.  The 
inverse  dependence  of  the  burning  rate  on  particle  size  implies  diffusion-controlled  combustion  ([122], 
[156]).  Equation  (57)  resembles  the  Mallard-Le  Chatelier  formula  [157]  for  the  flame  speed  of  a 
homogenous  gas-phase  mixture,  except  for  the  additional  term  accounting  for  the  energy  consumed  to 
vaporize  the  water. 

The  obtained  closed-form  expression,  Eq.  (57),  is  used  to  estimate  the  burning  rate  of  a 
stoichiometric  aluminum-water  mixture.  The  particle  size  and  pressure  are  taken  to  be  38  nm  and  3.65 
MPa,  respectively.  The  mean  specific  heat  of  the  mixture  is  taken  as  2.36  kJ/kg-K.  The  thermal 
conductivity  of  the  mixture  in  the  reaction  zone  is  calculated  to  be  1.21  W/m-K.  The  density  of  the 
mixture  is  assumed  to  be  equal  to  the  theoretical  value  of  1800  kg/m3.  The  enthalpies  of  reaction  and 
vaporization  of  water  are  taken  as  4400  and  1737  kJ/kg,  respectively.  The  burning  time  is  calculated  as 
0.07  ms.  Substituting  these  values  into  Eq.  (57),  the  burning  rate  is  estimated  to  be  7.28  cm/s,  which  is 
within  the  range  of  4.66-7.78  cm/s  observed  in  the  experiments  [138].  The  obtained  value  is  sensitive  to 
changes  in  the  burning  time  used  in  the  model.  For  example,  the  burning  rate  increases  from  7.28  to  8.6 
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cm/s,  when  the  burning  time  decreases  from  0.07  to  0.05  ms.  The  analysis  demonstrates  that  the  proposed 
model  properly  accounts  for  the  underlying  physicochemical  processes.  A  more  accurate  result  can  be 
obtained  by  relaxing  the  simplifying  assumptions. 


Table  7  Characteristics  of  aluminum  particles 


Particle  size,  nm 

Oxide  layer  thickness,  nm 

A1  content,  wt.% 

38 

3.10 

54.3 

50 

2.10 

68.0 

80 

2.70 

75.0 

120 

0.48 

96.5 

130 

2.20 

84.0 

Results  and  Discussion 

The  theoretical  framework  is  employed  to  calculate  the  flame  structure  and  burning  rate  of 
stoichiometric  mixtures  at  different  pressures  and  particle  sizes.  Table  7  shows  the  characteristics  of  the 
particles  considered  in  the  present  study.  The  thickness  of  the  oxide  layer  varies  in  the  range  of  0.48-3.10 
nm.  The  active  aluminum  content  decreases  with  decreasing  particle  size.  The  particle  composition 
significantly  influences  the  thermophysical  properties  of  the  mixture.  As  a  result,  it  is  important  to  use  an 
appropriate  value  of  the  oxide  layer  thickness  in  the  calculations.  Figure  48  shows  the  temperature 
distribution  for  a  stoichiometric  mixture  containing  38  nm  aluminum  particles  at  pressures  of  1  and  10 
MPa.  The  temperature  increases  from  an  initial  value  of  298  K  in  the  preheat  zone  and  attains  a 
maximum  value  of  1800  K  in  the  reaction  zone.  The  thickness  of  the  vapor  and  reaction  zones  decrease 
with  increasing  pressure.  The  temperature  distribution  is  further  altered  by  the  fact  that  the  vaporization 
temperature  of  water  increases  with  pressure.  An  estimate  of  the  reaction  zone  thickness  is  obtained  by 
multiplying  the  flame  propagation  velocity  and  particle  burning  time.  At  a  pressure  of  1  MPa,  the 
estimated  value  is  6  pm,  which  agrees  reasonably  well  with  the  result  of  the  present  analysis. 

Figures  (49)-(51)  show  the  effect  of  pressure  on  the  burning  rates  for  stoichiometric  mixtures 
containing  130,  80,  and  38  nm  particles,  respectively.  For  130  nm  particles,  the  burning  rate 


Figure  48  Effect  of  pressure  on  temperature  distribution  of  stoichiometric  A1-H20  mixture 
containing  38  nm  particles  in  (a)  the  preheat  zone;  (b)  the  reaction  zone 

increases  from  0.76  to  1.59  cm/s  when  the  pressure  increases  from  1  to  10  MPa.  A  similar  trend 
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is  observed  for  the  other  two  cases.  For  38  nm  particles,  the  presence  of  significant  scatter  in  the  burning 
rates  can  be  attributed  to  the  variations  in  the  packing  density  of  the  mixture.  The  actual  densities  are  in 
the  range  of  0.75-1.00  g/cm3,  which  are  lower  than  the  theoretical  value  of  1.80  g/cm3.  Such  a  disparity  is 
not  observed  for  80  and  130  nm  particles,  since  fewer  water  molecules  are  absorbed  on  the  particle 
surface.  The  burning  rate  decreases  with  increasing  packing  density  [138].  It  is,  thus,  not  surprising  that 
some  of  the  measured  burning  rates  are  greater  than  the  predicted  values.  The  pressure  exponent  in  the 
burning  rate  law  is  -0.3,  which  agrees  reasonably  well  with  experimental  data.  In  a  typical  composite 
solid  rocket  propellant,  the  burning  rate  increases  with  increasing  pressure.  A  general  explanation  for  this 
phenomenon  is  that  the  flame  stand-off  distance  decreases  with  increasing  pressure,  thereby  increasing 
the  heat  flow  to  the  propellant  surface  ([158],  [159]).  In  the  current  study,  a  visible  flame  appeared  to  be 
attached  to  the  burning  surface,  since  the  particles  undergo  heterogeneous  surface  reactions.  The  observed 
pressure  effect  can  be  attributed  to  the  fact  that  the  particle  burning  time  decreases  with 


Pressure,  MPa 

Figure  49  Effect  of  pressure  on  burning  rates  of  stoichiometric  A1-H20  mixtures  containing  130  nm 

particles. 


Figure  50  Effect  of  pressure  on  burning  rates  of  stoichiometric  A1-H20  mixtures  containing  80  nm 

particles 
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Pressure,  MPa 

Figure  51  Effect  of  pressure  on  burning  rates  of  stoichiometric  Al-  H20  mixture  containing  38  nm 

particles. 

with  pressure.  Figure  52  shows  the  effect  of  particle  size  on  the  burning  rate  for  a  stoichiometric  mixture 
at  a  pressure  of  3.65  MPa.  The  burning  rate  shows  a  particle  size  dependence  of  rb  =  a  dpn,  with  an 
exponent  of  -1.15.  The  experimental  data  suggest  that  the  burning  rate  is  inversely  proportional  to  the 
particle  diameter  [138].  In  the  present  analysis,  the  particle  burning  time  is  assumed  to  follow  dp2- law. 
The  actual  diameter  exponent  in  the  burning-time  relationship  may  be  slightly  lower  than  2,  when 
simultaneous  diffusion  of  water  vapor  and  aluminum  are  considered  [155].  This  may  explain  the  observed 
disparity  between  the  predicted  and  measured  values  of  the  diameter  exponent.  Figure  53  shows  a 
comparison  of  the  measured  and  calculated  burning  rates  with  those  obtained  using  the  following 
correlation: 


r6[cm/s]  =  98.8x(^[MPa])032(^[nm])-10  (58) 


Particle  Diameter  (dp),  nm 

Figure  52  Effect  of  particle  size  on  burning  rate  of  stoichiometric  A1-H20  mixture  at  3.65  MPa. 

Reasonably  good  agreement  is  achieved,  which  demonstrates  the  validity  of  the  correlation  for  situations 
encountered  in  the  current  study.  The  results  support  the  theory  that  the  rate-controlling  mechanism  is  the 
mass  diffusion  through  the  oxide  layers  of  the  particles.  In  the  present  analysis,  entrainment  of  particles  in 
the  gas  flow  has  been  neglected.  The  experiments  indicate  that  a  significant  number  of  particles  remained 
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in  the  quartz  tube  instead  of  being  transported  out  of  the  tube.  This  may  be  attributed  to  the  inertial  and 
gravitational  forces,  particle-particle  interactions,  and  quartz  tube  wall  effects.  The  region  in  which 
particle  motion  is  likely  to  be  important  is  the  post-combustion  zone,  since  the  particles  are  under  the 
continuous  influence  of  the  flow  of  the  combustion  gas  (H2).  This  can  be  incorporated  into  the  model  by 
considering  the  inertial  and  gravitational  forces,  interactions  and  collisions  between  particles,  and 
confining  effect  of  the  quartz  tube.  The  current  model  captures  the  main  features  of  alumimum-water 
combustion,  with  reasonably  good  agreement  between  experimental  data  and  model  predictions. 


Burning  Rate  (rb),  cm/s 

Figure  53  Measured  and  calculated  burning  rates  vs.  curve-fit  values  for  different  particle  sizes  and 
pressures,  r4[cm/s]  =  98.8 x (/?[MPa])° 3 V^nm])"1 0 

The  combustion  wave  propagation  of  nanoaluminum-water  mixtures  was  studied  theoretically 
and  experimentally  for  particles  in  the  size  range  of  38-130  nm  and  over  a  pressure  range  of  1-10  MPa. 
A  multi-zone  framework  was  established  to  predict  the  burning  and  flame  properties  by  solving  the 
conservation  equations  in  each  zone  and  imposing  mass  and  energy  continuities  at  the  interfacial 
boundaries.  The  flame  propagation  characteristics  were  measured  by  burning  nanoaluminum-water 
strands  in  a  constant-volume  vessel.  Emphasis  was  placed  on  the  effects  of  particle  size  and  pressure. 
An  analytical  expression  for  the  burning  rate  was  derived,  and  physicochemical  parameters  that  dictate 
the  flame  behavior  were  identified.  For  conditions  present  in  the  study,  the  burning  rate  showed  pressure 

and  particle  size  dependencies  of  ^[cm/s]  =  98.8x(^[MPa])032(^[nm])_1 0 .  The  flame  thickness 

increased  with  increasing  particle  size  and  decreasing  pressure.  Results  supported  the  hypothesis  that  the 
combustion  of  aluminum-water  mixtures  is  controlled  by  mass  diffusion  across  the  oxide  layers  of  the 
particles.  Future  work  will  specifically  examine  the  effects  of  particle  entrainment  and  agglomeration  on 
burning  properties  of  nanoaluminum-water  mixtures. 

C.  Pyrophoricity  of  Nano-Aluminum  Particles 

Aluminum  particles  are  extensively  used  in  propulsion  and  energy-conversion  applications  due  to 
their  high  oxidation  enthalpy,  relative  safety,  and  low  cost  [160].  When  a  nascent  aluminum  particle  is 
exposed  to  an  oxidizing  environment  at  room  temperature,  the  chemical  reactions  between  the  aluminum 
atoms  and  oxidizer  molecules  result  in  the  formation  of  an  amorphous  aluminum  oxide  (A1203)  layer  on 
the  particle  surface  [161].  Typically,  the  rate  of  heat-release  due  to  these  chemical  reactions  is  lower  than 
that  of  heat  loss  to  the  surrounding  environment.  As  a  result,  the  particle  does  not  ignite  and  the  thickness 
of  the  oxide  layer  saturates  at  a  value  in  the  range  of  0.5-4  nm,  depending  on  the  temperature  and  duration 
of  exposure  to  the  oxidizing  gas  [161].  The  oxide  layer  protects  the  aluminum  core  from  further  attack  of 
the  oxidizer  molecules.  The  ignition  temperature  of  passivated  aluminum  particles  decreases  from  2350  K 
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(melting  point  of  the  oxide  layer)  at  100  pm  to  933  K  (melting  point  of  the  aluminum  core)  at  0.1  pm 
[162].  The  enhanced  reactivity  of  nano-sized  particles  poses  significant  safety  hazards.  If  the  size  of  the 
particle  is  decreased  below  a  critical  value,  the  particle  could  ignite  when  exposed  to  an  oxidizing  gas  at 
room  temperature,  a  phenomenon  known  as  pyrophoricity  [163].  It  is  a  major  safety  issue  during  particle 
manufacture,  handling,  and  storage.  It  can  also  be  employed  for  useful  applications  like  decoy  flare  for 
defending  heat-seeking  missiles  [164].  Reliable  measurements  or  predictions  of  the  critical  particle  size 
for  pyrophoricity  are  of  paramount  interest. 

Theoretical  studies  on  metal  pyrophoricity  are  very  limited.  Glassman  et  al.  [163]  postulated  that 
a  metal  particle  is  pyrophoric  if  the  chemical  energy  release  is  sufficient  to  vaporize  the  metal.  A  quasi 
steady-state  energy  balance  was  performed.  For  an  oxide-layer  thickness  of  2.5  nm,  the  critical  particle 
diameter  was  calculated  to  be  23  nm.  The  study  ignored  the  effect  of  heat  losses  to  the  ambient 
environment  and  considered  bulk  material  properties.  For  nano-sized  particles,  the  melting  and  boiling 
points  [165-171],  heat  of  fusion  [168-171],  vapor  pressure  [172-174],  ignition  temperature  [3],  and  heat  of 
reaction  [175]  are  substantially  different  from  their  respective  bulk  values.  Puri  [176]  estimated  the 
critical  particle  size  by  performing  transient  energy  balance  in  conjunction  with  the  use  of  size-dependent 
thermophysical  properties.  The  oxidation  was  assumed  to  be  controlled  by  diffusion  of  oxygen  molecules 
through  the  gas-phase  mixture  to  the  particle  surface.  The  estimated  critical  particle  size  is  20  nm.  It  is 
worth  noting  that  once  a  monomolecular  oxide  layer  is  formed,  the  oxidation  of  nano-aluminum  particles 
is  governed  by  the  Mott-Cabrera  oxidation  kinetics  [177,178].  An  essential  feature  of  this  model  is  that 
the  metal  electrons  transverse  the  thin  oxide  layer  either  by  thermionic  emission  or  tunneling.  The 
electrons  ionize  the  adsorbed  oxygen  atoms  to  create  an  electrostatic  potential  between  the  oxide-oxidizer 
and  oxide-metal  interfaces.  Figure  54  shows  the  induced  electric  field  in  the  particle,  which  significantly 
lowers  the  energy  barrier  for  metal  ion  diffusion,  resulting  in  higher  oxidation  rates.  Recently,  the  Mott- 
Cabrera  oxidation  kinetics  was  employed  to  study  the  pyrophoricity  of  nascent  aluminum  particles  [179]. 
The  critical  particle  diameter  was  estimated  to  be  68  nm,  which  is  significantly  greater  than  the  values 
predicted  by  Glassman  et  al.  [163]  and  Puri  [176].  The  study  did  not  consider  the  size  dependence  of 
physicochemical  properties.  The  sensitivity  of  the  results  to  changes  in  the  model  parameters  has  not  been 
investigated.  In  addition,  the  effect  of  the  passivating  oxide  layer  on  pyrophoricity  of  particles  needs  to 
be  examined. 


oxidizer 


Figure  54  Electric-field  across  the  oxide  layer  in  a  spherical  nano-particle  surrounded  by  oxidizer 
molecules  (+  positive  charge,  -  negative  charge)  according  to  the  Mott-Cabrera  mechanism  [177]. 

In  the  present  study,  a  more  comprehensive  analysis  is  performed  by  taking  into  account  transient 
energy  balance  along  with  accurate  evaluation  of  physicochemical  properties  of  nano-sized  particles.  The 
oxidation  mechanism  is  based  on  a  recently  developed  Mott-Cabrera  oxidation  kinetics,  which  is 
applicable  for  spherical  nano-particles  [178].  The  sensitivity  of  the  model  results  to  the  choice  of 
physicochemical  properties  of  the  particle  and  gas,  polymorphic  state  of  the  oxide  layer,  parameters  of  the 
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Mott-Cabrera  oxidation  kinetics,  and  heat-transfer  correlation  is  analyzed.  Both  nascent  and  passivated 
particles  are  considered.  The  temperature  and  pressure  of  the  ambient  environment  are  chosen  as  300  K 
and  1  atm,  respectively.  Results  from  the  present  work  are  compared  with  those  of  previous  theoretical 
and  experimental  studies.  Reasonably  good  agreement  with  experimental  data  is  achieved,  thereby 
demonstrating  the  validity  of  the  model. 


Physicochemical  Properties  of  Nano-Sized  Aluminum  Particles 


The  physicochemical  properties  of  nano-sized  aluminum  particles  are  significantly  different  from 
their  counterparts  of  the  bulk  material  [21].  The  fraction  of  atoms  on  the  surface  layer  of  an  aluminum 
particle  is  given  by 
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where  8  =  286  pm  and  D  the  particle  diameter.  It  increases  from  2  to  92  %,  when  the  particle  size 
decreases  from  100  to  1  nm.  The  surface  atoms  have  lower  coordination  numbers  and  higher  energies 
than  the  atoms  in  the  interior  of  the  particle,  thereby  leading  to  size-dependent  properties.  For  example, 
the  melting  temperature  of  a  nascent  aluminum  particle  decreases  from  937  K  at  10  nm  to  473  K  at  2  nm 
[165].  The  melting  point,  Tm,  and  the  enthalpy  of  fusion,  Lfus ,  of  aluminum  particles  can  be  written  as 
[171] 
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where  R  is  the  particle  radius,  Ct  =  1.9186,  and  R0  =  0.9492  nm.  The  subscript  b  refers  to  the  bulk 
material. 

Another  thermophysical  property  of  relevance  to  the  present  study  is  the  boiling  temperature, 
defined  as  the  temperature  at  which  the  vapor  pressure  is  equal  to  the  ambient  pressure.  The  vapor 
pressure  of  aluminum  near  the  surface  of  a  liquid  droplet,  pD,  is  calculated  using  the  Kelvin  equation 
[172] 

Pd  =  Po  exp(4crW  kBTD).  (62) 


where  D  is  the  diameter  of  the  particle,  k!;  the  Boltzmann  constant,  and  v  the  molar  volume.  The  vapor 
pressure  of  aluminum  over  a  flat  surface,  and  the  surface  tension  of  aluminum,  a,  are  calculated 
respectively  as  follows  [173,174]: 
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cr  =  948-0.2027, 


(64) 


where  p  is  the  pressure  and  T  the  temperature.  Equations  (62)-(64)  are  solved  iteratively  to  calculate  the 
boiling  temperature  as  a  function  of  the  particle  size.  Figure  55  shows  the  result  at 
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Figure  55  Vaporization  temperature  of  aluminum  as  a  function  of  particle  size  over  a  pressure 

range  of  0.5-2.0  atm. 


three  different  pressures  of  0.5,  1,  and  2  atm.  The  boiling  temperature  of  aluminum  increases  with 
increasing  pressure.  It  decreases  significantly  from  the  bulk  value  for  particles  smaller  than  10  nm,  a  trend 
that  is  consistent  with  those  observed  for  the  melting  temperature  and  enthalpy  of  fusion. 

The  heat  of  reaction  is  another  physicochemical  property  that  needs  to  be  specified.  It  determines 
the  energy  release  during  the  formation  of  the  oxide  layer  on  a  nascent  particle.  The  theoretical  data 
available  in  Ref.  [175]  is  curve-fitted  to  obtain  a  correlation  for  the  heat  of  reaction  of  aluminum  particles 
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where  AHr  is  the  heat  of  reaction  and  D  the  particle  size  in  nm.  The  heat  of  reaction  decreases  with 
decreasing  particle  size.  In  other  words,  the  oxidation  reaction  is  less  exothermic  for  smaller  particles. 
This  effect  becomes  more  substantial  for  particles  smaller  than  10  nm. 

Another  parameter  of  concern  is  the  density  of  the  aluminum  oxide,  which  exists  in  many 
metastable  polymorphs.  The  hexagonally  packed  a-Al203  (corundum)  has  a  density  of  4  g/cm3  [161]. 
Other  polymorphs  such  as  y,  S. [  and  # phases  have  face-centered  cubic  (FCC)  structures.  The  initial  oxide 
layer  covering  the  aluminum  particle  is  amorphous  [161].  The  density  of  the  amorphous  oxide  layer  is 
taken  to  be  3.05  g/cm3  [161]. 

Steady  Homogenous  Gas-Phase  Reaction  Theory 

The  analysis  presented  in  Ref.  [163]  assumes  homogeneous  gas-phase  combustion  of  aluminum 
vapor  and  oxidizing  species.  It  neglects  heat  losses  to  the  ambient  environment  and  size-dependence  of 
physicochemical  properties.  Under  these  assumptions,  the  critical  condition  for  metal  pyrophoricity  was 
obtained  by  equating  the  chemical  heat  release  to  the  sum  of  the  energies  needed  to  heat  the  particle  to  its 
boiling  point  and  vaporize  the  metal  [163]: 
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Particle  Diameter,  nm 

Figure  56  Energy  distribution  for  the  oxidation  of  aluminum  particles  with  an  oxide  layer  thickness 
of  2.5  nm:  steady  homogeneous  gas-phase  reaction  theory  with  bulk  material  properties 

where  S  is  the  oxide  layer  thickness,  p  the  density,  Tb  the  boiling  point  of  aluminum,  Lvap  the  enthalpy  of 
vaporization  of  aluminum,  and  H  the  enthalpy.  The  subscripts  ox  and  cr  refer  to  oxidation  and  critical 
condition,  respectively.  For  the  sake  of  comparison,  bulk  values  of  thermophysical  properties  are 
employed.  They  are  obtained  from  the  JANNAF  tables  and  are  given  in  Table  8  [181].  Figure  56  shows 
the  calculated  energy  budget  of  the  oxidation  process  for  different  particle  sizes.  The  thickness  of  the 
oxide  layer  is  chosen  as  2.5  nm.  The  curves  corresponding  to  the  total  absorbed  and  released  energies 
intersect  at  a  particle  size  of  28  nm.  The  predicted  critical  particle  size  agrees  reasonably  with  the  value  of 
28.6  nm  reported  in  Ref.  [163],  thereby  demonstrating  the  accuracy  of  the  present  calculation. 

The  size-dependent  physicochemical  properties  are  used  to  obtain  a  revised  estimate  of  the 
critical  particle  size.  The  density  of  the  oxide  layer  is  taken  as  3.05  g/cm3 

Table  8  Thermo-physical  properties  of  bulk  aluminum  and  alumina 


Property 

Value 

Density  of  alumina 

Density  of  aluminum 

Heat  release  during  oxidation,  ~H29^ox 

4000  kg/m3 
2700  kg/m3 
1675  kJ/mol 

TT°  IT o 

21  Tb  21  298 

454  kJ/mol 

H°Tb-H°m+Lv 

381.67  kJ/mol 

Aluminum  melting  point 

940  K 

Aluminum  boiling  point 

2740  K 

Alumina  melting  point 

2327  K 

Alumina  volatilization  temperature 

4000  K 
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Oxide  Layer  Thickness  (8),  nm 

Figure  57  Critical  pyrophoricity  ratio  and  core  diameter  as  a  function  of  the  oxide  layer  thickness: 
steady  homogeneous  gas-phase  reaction  theory  with  size-dependent  properties 

[161].  Figure  57  shows  the  variation  of  the  critical  pyrophoricity  ratio,  {Sj 7?)  ^  and  particle  size  with  the 

oxide  layer  thickness.  The  former  decreases  with  increasing  oxide  layer  thickness.  For  an  oxide  layer 
thickness  of  2.5  nm,  it  takes  a  value  of  0.217.  The  corresponding  critical  particle  diameter  (inclusive  of 
the  oxide  layer)  is  23  nm,  which  is  lower  than  the  value  of  28  nm  predicted  using  the  constant-property 
assumption.  The  effect  of  pressure  on  the  critical  particle  size  is  also  studied  in  the  range  of  0.5-2  atm 
Pressure  plays  a  negligible  role  in  the  range  of  0.5-2  atm.  Note  that  the  above  estimates  only  serve  as  a 
guideline,  since  the  analysis  neglects  the  effects  of  heat  losses  and  kinetics  of  oxidation  of  aluminum 
particles 

Transient  Heterogeneous  Surface  Reaction  Theory 

When  a  nascent  particle  is  exposed  to  the  oxidizing  gas,  the  temperature,  oxide  layer  thickness, 
and  core  diameter  vary  continuously  with  time.  As  a  result,  the  steady-state  assumption  is  not  valid. 
Particle  ignition  is  a  transient  process,  and  the  corresponding  ignition  delay  can  be  characterized  by  the 
point  at  which  temperature  runaway  occurs.  During  ignition,  the  rate  of  chemical  energy  release  competes 
with  that  of  heat  loss  to  the  oxidizing  gas.  An  unsteady  energy  balance  is  thus  required  to  obtain  an 
accurate  estimate  of  the  critical  particle  size  for  pyrophoricity.  The  analysis  follows  the  approach  adopted 
by  Mohan  et  al.  [179],  but  is  extended  to  include  the  size  and  temperature  dependencies  of  particle 
physicochemical  properties.  Application  of  the  energy  conservation  leads  to  the  following  equation. 

^  dT  , 

mCp~dt=h°x~q°un  (67) 

where  m  is  the  mass  of  the  particle,  Cp  the  specific  heat  of  the  particle,  T  the  temperature,  t  the  time,  h()x 
the  rate  of  chemical  energy  release.  The  density  and  specific  heat  are  calculated  as  a  function  of  the 
particle  temperature.  The  rate  of  heat  loss  to  the  surrounding  gas,  qout ,  consists  of  contributions  from 
conduction  and  radiation,  denoted  by  the  subscripts  c  and  r,  respectively,  as  follows: 
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Qou.  =4c+4r, 

In  the  continuum  regime,  the  conductive  heat  transfer  between  the  particle  and  surrounding  gas  is 
controlled  by  energy  diffusion,  given  by  [182] 

qc=4xR\(T-Te).  (69) 

where  R  is  the  particle  radius  and  Xg  the  thermal  conductivity  of  the  gas.  The  subscript  e  refers  to  the 
ambient  condition.  At  a  temperature  of  1000  K  and  pressure  of  1  atm,  the  mean  free  path  of  oxygen 
molecules  is  223  nm,  which  is  significantly  greater  than  the  particle  size.  As  a  result,  the  continuum  heat 
transfer  correlation  ceases  to  be  valid.  In  the  free-molecular  regime,  the  conductive  heat-flux  is  replaced 
by  the  molecular  heat  flux,  qmol ,  given  by  [183]: 
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where  M  is  the  mass  of  the  oxygen  molecule,  7*  the  adiabatic  constant  calculated  at  a  temperature  T*  = 
(T+Te)/2,  and  a  =  0.85  the  accommodation  coefficient  [184].  The  radiation  heat  transfer  can  be  calculated 
using  the  Stefan-Boltzmann  relation: 


qr  =sA<j[Ta  (7i) 

where  A  is  the  particle  surface  area,  a  the  Stefan-Boltzmann  constant,  and  s  =  0.2  the  emissivity 
of  the  oxidized  aluminum  surface  [185].  Equation  (13)  is  strictly  not  valid  for  nano-sized  particles. 


Particle  Diameter,  nm 

Figure  58  Effect  of  particle  size  on  emissivity  of  aluminum  nano-particles  at  different  temperatures 
obtained  using  the  magnetic  dipole  approximation  theory  [186] 


The  size  and  temperature  dependencies  of  the  particle  emissivity  must  be  considered.  Figure  58  shows 
the  effect  of  particle  size  on  the  emissivity  of  aluminum  particles  at  different  temperatures  calculated 
using  the  magnetic  dipole  approximation  theory  given  in  Ref.  [186].  The  emissivity  of  nano-sized 
aluminum  is  several  orders  of  magnitude  lower  than  its  bulk  value.  A  more  accurate  value  of  the  particle 
emissivity  can  be  obtained  by  considering  the  electrical  dipole  contribution,  which  may  become 
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significant  at  higher  temperatures  [187].  The  analysis,  however,  indicates  that  the  model  results  are 
insensitive  to  changes  in  the  particle  emissivity. 

The  heat-release  term  is  calculated  using  the  Mott-Cabrera  oxidation  kinetics  [188].  Ermoline  and 
Dreizin  [189]  incorporated  the  electric-field  correction  and  effects  of  volume  changes  in  a  growing  oxide 
layer  and  a  shrinking  aluminum  core,  and  derived  the  following  equations  for  spherical  particles: 

Table  9  Constants  in  Mott-Cabrera  equations 


Constant  Value 


n 

10  nm'2 

V 

1012  s'1 

a 

0.12  nm 

( Pm, 

-1.6  V 

w 

2.6  eV 

q 

3e 

-0.0166  nnr 

Q2 

0.023  nm3 

dS_ 

dt 

dK 

dt 


(Qj  +Q2) 
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kRj 
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V  kBT  j 
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exp 
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where  Rc  is  the  core  radius  and  8  the  oxide  layer  thickness.  The  model  constants  are  given  in  Table  9. 
The  rate  of  chemical  heat  release  can  be  written  as 


t"  A  9  T  tiiV 

hox=^Rc  PaiK-t-’  (73) 

dt 

where  hr  is  the  heat  of  reaction.  Equations  (72)  and  (73)  require  the  presence  of  an  oxide  layer  on  the 
particle.  The  formation  of  the  monomolecular  oxide  layer  on  a  nascent  particle  is  extremely  fast  and  can 
be  treated  as  an  adiabatic  process  [190].  A  detailed  justification  is  given  in  Ref.  [190].  Our  calculations 
indicate  that  the  growth  of  the  oxide  layer  is  adiabatic  up  to  a  thickness  of  0.3  nm.  This  further 
demonstrates  the  validity  of  the  assumption.  As  a  result,  chemical  equilibrium  analysis  can  be  performed 
to  calculate  the  particle  temperature  upon  the  formation  of  the  0.3-nm  thick  (monomolecular)  oxide  layer. 
No  such  calculation  is  necessary  for  passivated  particles.  The  energy  balance  based  on  chemical 
equilibrium  analysis  can  be  written  as 

(mAiCp,Al+moXCp,ox)(Tf  -300 )  +  mAlLfusAl=-moxhox,  (74) 

where  7}  is  the  final  temperature. 

Results  and  Discussion 

The  theoretical  framework  described  in  Section  4  is  used  to  analyze  the  pyrophoricity  of  nascent 
and  passivated  aluminum  particles  at  nano  scales.  Initially,  the  nascent  particle  is  assumed  to  be  devoid  of 
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the  oxide  layer.  For  passivated  particles,  the  initial  values  of  the  oxide  layer  thickness  in  the  range  of  0.3- 
0.5  nm  are  considered. 

Nascent  aluminum  particles 

The  growth  of  an  oxide  layer  on  a  nascent  aluminum  particle  is  divided  into  two  stages.  In  the 
first  stage,  a  0.3-nm  thick  monomolecular  oxide  layer  is  formed.  A  chemical  equilibrium  analysis  is 
performed  to  calculate  the  particle  temperature  after  the  formation  of  the  oxide  monolayer.  A  companion 
transient  energy  balance  analysis  employing  the  Mott-Cabrera  kinetics  is  conducted  to  determine  the 
minimum  ignition  temperature  for  particles 
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Figure  59  Variations  of  particle  temperature  and  oxide  thickness  with  time  for  a  core  diameter  of 

10  nm  and  oxide  layer  thickness  of  0.3  nm. 

encapsulated  with  a  0.3-nm  thick  oxide  layer.  Ignition  is  assumed  to  be  achieved  when  the  particle 
temperature  increases  monotonically  beyond  the  melting  point  of  the  oxide  film.  If  the  ignition 
temperature  is  lower  than  the  result  of  chemical  equilibrium  analysis,  the  particle  can  self-ignite  due  to 
heterogeneous  oxidation  reactions  beginning  at  300  K.  It  is,  thus,  considered  to  be  pyrophoric.  Figure  59 
shows  the  variations  of  the  temperature  and  oxide  layer  thickness  with  time  for  a  particle  with  a  core 
diameter  of  10  nm.  For  an  initial  temperature  of  680  K,  the  rate  of  chemical  heat-release  is  not  sufficient 
to  balance  the  heat  loss  to  the  ambient  environment.  As  a  result,  temperature  runaway  does  not  occur.  The 
isotherm  corresponds  to  the  melting  of  the  core  at  750  K,  which  is  lower  than  the  bulk  melting  point  of 
aluminum  by  183  K.  The  oxide  layer  thickness  follows  a  similar  trend.  When  the  initial  temperature  is 
chosen  as  690  K,  the  temperature  and  oxide  layer  thickness  increase  sharply  after  a  time  period  of  20  ns. 
The  rate  of  chemical  energy  release  is  significantly  greater  than  that  of  heat  loss  to  the  ambient 
environment.  The  same  phenomenon  takes  place  for  all  temperatures  greater  than  690  K.  For  a  particle 
with  a  10  nm  core,  the  minimum  ignition  temperature  is,  thus,  taken  as  690  K.  The  sensitivity  of  the 
results  to  the  ambient  temperature  stems  from  the  fact  that  the  reaction  rate  bears  an  exponential 
dependence  on  temperature.  A  similar  analysis  is  performed  for  different  particle  sizes  and  the 
corresponding  minimum  ignition  temperature  is  determined. 

Figure  60  shows  the  comparison  of  the  results  of  chemical  equilibrium  and  transient  energy 
balance  analyses  for  two  different  cases.  In  the  first  case,  bulk  material  properties  calculated  at  room 
temperature  are  employed.  A  crystalline  aluminum  oxide  layer  with  a  density  of  4.0  g/cm3  is  assumed  to 
cover  the  active  aluminum  surface.  In  the  second  case,  size  and  temperature  dependent  properties  are  used 
and  an  amorphous  oxide  layer  with  a  density  of  3.05  g/cm3  is  considered.  The  chemical  equilibrium 
analysis  indicates  that  the  final  particle  temperature  decreases  with  increasing  particle  size,  since  more 
energy  is  spent  to  heat  a  larger  particle.  The  transient  energy  balance  analysis,  on  the  other  hand,  suggests 
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that  the  minimum  ignition  temperature  increases  with  increasing  particle  size,  since  the  rate  of  heat  loss  to 
the  oxidizing  gas  is  proportional  to  the  particle  surface  area.  The  plateau-like  feature  indicates  that  the 
corresponding  particles  need  to  be  pre-heated  to  the  melting  point  of  the  core.  These  two  curves  intersect 
at  the  critical  particle  size,  the  specific  value  of  which  is  different  for  the  two  cases.  For 


Figure  60  Particle  temperature  upon  growth  of  0.3  nm  oxide  layer  and  minimum  temperature 
necessary  for  ignition  for  (a)  crystalline  oxide  layer  and  bulk  material  properties  at  room 
temperature;  (b)  amorphous  oxide  layer  with  temperature  and  size  dependent  material 

the  first  case,  the  calculated  value  is  45  nm,  which  is  comparable  with  the  estimate  reported  in  Ref.  [190]. 
For  the  second  case,  the  critical  particle  size  is  predicted  to  be  32  nm.  The  model  results  were  sensitive  to 
the  changes  in  the  parameters  of  the  Mott-Cabrera  oxidation  kinetics.  For  example,  increasing  the  Mott 
potential  from  1.6  to  1.75  V  resulted  in  an  increase  in  the  critical  particle  size  by  18  %.  In  the  above 
analyses,  the  free-molecular  regime  heat  transfer  model  is  employed.  Figure  61  shows  the  result  when  the 
continuum  hypothesis  is  invoked.  It  overestimates  the  heat  losses  to  the  ambient  environment.  As  a  result, 
particles  need  to  be  heated  to  higher  temperatures  for  ignition  to  occur.  The  calculated  critical  particle 
size  of  18  nm  is  significantly  lower  than  the  value  obtained  using  the  model  for  the  free-molecular 
regime.  The  continuum  heat  transfer  correlation  gives  inadequate  description  of  ignition  and  combustion 
of  nano-sized  particles.  Table  10  shows  the  comparison  between  the  predictions  of  different  models 


Figure  61  Critical  particle  size  predicted  using  the  continuum  heat  transfer  correlation 
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and  experimental  data.  Results  from  the  present  analysis  exhibit  reasonably  good  agreement  with 
measurements. 


Table  10  Critical  particle  size  predicted  by  different  models  and  their  comparison  with  the 

experimental  data 


Model 

Critical  particle  size,  nm 

Steady-state, 
homogeneous  gas- 

phase  reaction  theory 

Constant  properties  [163] 

23 

Size-dependent  properties 

18 

Transient 

heterogeneous  surface 
reaction  theory 

Diffusion  oxidation  model  [176] 

20 

Mott-Cabrera 

kinetics 

Free  molecular 

32 

Continuum 

18 

Experimental  data  [163] 

30 

In  the  above  analysis,  thermal  accommodation  coefficient  was  taken  as  0.85.  In  reality,  the 
accommodation  coefficient  is  a  temperature  dependent  parameter  [187].  Figure  62  shows  the  result  when 
temperature  dependence  of  the  thermal  accommodation  coefficient  is  considered.  The  accommodation 
coefficient  decreases  with  increasing  temperature  and  a  result  the  critical  particle  size  increases  from  32 
nm  to  40  nm. 


Figure  62  Critical  particle  size  predicted  using  the  free-molecular  heat  transfer  correlation  in 
conjunction  with  temperature  dependent  accommodation  coefficient. 

Passivated  aluminum  particles 

The  same  approach  outlined  previously  is  adopted  for  passivated  aluminum  particles,  except  that 
the  chemical  equilibrium  analysis  is  not  employed.  The  oxide  layer  thickness  is  varied  in  the  range  of  0.3- 
0.5  nm.  Figure  63  shows  the  temporal  evolution  of  the  temperature  for  particles  with  a  core  diameter  of 
3.8  and  6  nm  along  with  an  oxide  layer  thickness  of  0.3  nm.  At  300  K,  temperature  runaway  is  observed 
only  when  the  particle  size  is  3.8  nm  and,  thus,  the  critical  particle  size  is  3.8  nm.  Figure  64  shows  the 
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temperature  evolution  for  a  core  diameter  of  3.8  nm  and  shell  thickness  in  the  range  of  0.3-0. 5  nm. 
Particles  with  0.4  and  0.5  nm  thick  oxide  layers  remain  stable  at  300  K,  suggesting  that  oxide  layers 
thicker  than  0.3  nm  render  the  particle  to  be  non-pyrophoric.  Aluminum  particles  are  typically  covered  by 
an  oxide  layer,  which  can  have 


Time,  ns 

Figure  63  Variation  of  particle  temperature  with  time  for  core  diameters  of  3.8  and  6  nm  and  oxide 

layer  thickness  of  0.3  nm. 

thickness  in  the  range  of  0.5-4  nm  [161].  Moreover,  particles  smaller  than  3.8  nm  are  not  of  practical 
interest.  As  a  result,  commercially  available  passivated  particles  are  found  to  be  non-pyrophoric. 


Figure  64  Variation  of  particle  temperature  with  time  for  a  core  diameter  of  3.8  nm  and  oxide  layer 

thickness  of  0.3,  0.4,  and  0.5  nm 


Pyrophoricity  of  nano-sized  nascent  and  passivated  aluminum  particles  in  air  was  studied  based 
on  a  transient  energy  balance  analysis.  The  work  considered  radiation  and  molecular  heat  exchange  with 
the  surrounding  environment  and  employed  the  Mott-Cabrera  oxidation  mechanism  for  the  heterogeneous 
chemical  reactions  that  occur  on  the  particle  surface.  The  size-dependence  of  physicochemical  properties 
of  the  particle  is  also  incorporated.  Results  suggested  that  nascent  aluminum  particles  with  diameters 
lower  than  32  nm  are  pyrophoric.  The  critical  particle  size  for  particles  passivated  with  a  0.3  nm  thick 
oxide  layer  was  calculated  as  3.8  nm.  Commercially  available  nano-scale  aluminum  particles  with  an 
oxide-layer  thickness  in  the  range  of  0. 5-4.0  nm  were  predicted  to  be  non-pyrophoric.  The  model  results 


75 


were  sensitive  to  the  choice  of  physicochemical  properties,  the  Mott-Cabrera  oxidation  parameters,  and 
the  polymorphic  state  of  the  oxide  layer.  The  continuum  hypothesis  gave  inadequate  description  of 
pyrophoricity  of  aluminum  particles  at  nano  scales. 

IV.  Experimental  Studies  on  Nanometer  Aluminum  and  Composite  Micron  Particles 
with  Nanometer  Features  (Purdue  University) 

A.  Studies  of  the  reaction  of  Aluminum  with  H20  and  H202 

Over  the  past  few  years,  the  combustion  of  nano-aluminum/water  (nAl/  H20)  propellants  has 
been  widely  reported,  but  further  progress  has  been  slowed  for  the  following  reason:  the  loosely 
correlated  trends  in  combustion  data  are  insufficient  in  guiding  further  research  efforts,  and  they  cannot 
be  used  to  significantly  improve  the  Isp  observed  in  static  rocket  motor  tests.  It  was  previously  found  that 
different  mixing  techniques  (hand,  planetary  and  resonant  mixers,  duration  and  temperature),  or 
equivalence  ratio  gave  rise  to  different  burning  rates,  but  the  influence  of  pH  and  rheology  on  nAl/  H20 
propellants  was  not  considered.  We  have  found  in  our  studies  that  the  effects  of  pH  on  nAl/  H20 
propellants  are  profound,  and  correlate  well  with  viscosity,  low  pressure  deflagration  limits,  burning  rate 
exponents,  and  rocket  motor  performance.  Our  findings  suggested  that  coagulation  can  influence  the 
pressure  exponent  over  a  wide  range  of  values  (0.34-0.68).  For  particle  diameters  <1  pm,  dispersion 
during  mixing  is  affected  more  by  electrostatic  repulsion  from  charged  ions  than  from  mechanical 
agitation,  and  this  is  reflected  through  zeta  potential  and  viscosity  measurements  at  different  pH  levels. 
Additionally,  we  observed  that  pH  has  an  influence  on  nAl/H20  reaction  kinetics  during  ignition,  as  the 
propellant  transitions  from  low  temperature  oxidation  to  high  temperature  combustion.  See  DOI: 
10.1080/00102202.2012.759948  for  more  details.  Also  see  DOI:  10.1002/prep.201200143  for  other  work 
characterizing  safety,  aging  and  equivalence  ratio  effects  on  nAl  and  water  reactions. 

We  also  completed  studied  involving  the  burning  rate  characteristics  of  hydrogen  peroxide  and 
micron-aluminum  propellants.  Theoretical  calculations  show  that  the  sea  level  specific  impulse  of  this 
simple  binary  mixture  is  comparable  to  standard  composite  propellant.  In  addition,  the  aluminum  particle 
size,  hydrogen  peroxide  concentration,  and  mixture  ratio  can  be  adjusted  over  a  flat  peak  performance 
regime  to  attain  specific  thrust  profiles  and  durations.  We  measured  the  burning  rates  in  a  windowed 
pressure  vessel  at  pressures  ranging  from  7  to  14  MPa.  Results  show  that  mixture  burning  rates  span  from 
0.5  to  4.5  cm/s  at  7  MPa  with  power  law  burning  rate  pressure  exponents  ranging  from  0.33  to  1.07.  In 
this  study,  we  focus  a  statistical  analysis  on  the  determination  of  the  most  influential  variables  affecting 
the  burning  rate  and  applied  a  thermal  analysis  to  determine  the  combustion  regimes  of  these  mixtures. 
The  statistical  analysis  provided  a  multivariate  regression  model  for  the  logarithm  of  the  burning  rate  with 
a  correlation  coefficient  of  0.93.  The  model  suggests  aluminum  diameter  is  the  most  important  factor 
affecting  the  overall  burning  rate,  and  H202  concentration  as  the  most  influential  variable  on  the  burning 
rate  pressure  dependence.  The  burning  rate  dependence  on  theoretical  combustion  temperature  shows  two 
distinct  combustion  regimes  attributed  to  kinetic  and  diffusion  controlled  combustion.  A  simple  thermal 
analysis  confirms  the  experimental  burning  rate  pressure  dependence  observed  for  these  two  regimes.  See 
DOI:  10. 1016/j.combustflame. 2012. 10.001  for  more  details. 

B.  Encapsulating  Nanoscale  Particles  in  Oxidizer  Crystals 

The  addition  of  nanoparticles  to  composite  propellant  formulations  theoretically  leads  to 
dramatically  increased  burning  rate  and  combustion  efficiency.  However,  nanoparticles  have  played  a 
severely  limited  role  in  modem  propellant  formulation  development  thus  far  because  of  undesirable 
resulting  properties.  The  literature  is  replete  with  examples  of  unworkable  viscosities  during  processing, 
brittle  finished  products,  poor  combustion  efficiency,  and/or  combustion  stability  problems  during  use. 
These  issues  ultimately  are  a  result  of  the  high  surface  area  (25-250  m2g-l)  of  the  nanoparticles 
preferentially  adsorbing  the  liquid  propellant  binder  components  during  processing,  preventing  complete 
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dispersion  and  reducing  the  capability  of  the  binder  matrix  to  retain  other  propellant  ingredients. 
Including  nanoscale  Energetic  crystals  with  nanoparticle  inclusions  were  prepared  and  characterized, 
focusing  on  the  nanosized  iron(III)  oxide-ammonium  perchlorate  system,  generated  using  an  ethyl 
acetate/acetone  antisolvent:  solvent  system.  We  have  shown  that  capture  is  dependent  on  antisolvent-to- 
solvent  ratio;  increased  quantities  of  antisolvent  lead  to  faster  nucleation  rates,  smaller  crystals,  and 
improved  capture.  Additionally,  the  crystal  habit  formation  is  modified  by  the  addition  of  nanoparticles, 
resulting  in  highly  uniform  rod-like  crystal  structures.  Results  with  nanoaluminum-ammonium 
perchlorate  and  nanoaluminum-cyclotrimethylenetrinitramine  (RDX)  systems  have  shown  similar  capture 
behavior.  See  DOI:  10.1002/prep.201200142  for  more  details. 

C.  Inclusion  Materials  in  Micron  Scale  Aluminum 

Micrometer-sized  aluminum  is  widely  used  in  energetics;  however,  performance  of  propellants, 
explosives,  and  pyrotechnics  could  be  significantly  improved  if  its  ignition  barriers  could  be  disrupted. 
We  report  morphological,  thermal,  and  chemical  characterization  of  fuel  rich  aluminum 
polytetrafluoroethylene  (70-30  wt-%)  reactive  particles  formed  by  high  and  low  energy  milling.  Average 
particle  sizes  range  from  15-78  pm;  however,  specific  surface  areas  range  from  approx.  2-7  m2g-l  due  to 
milling  induced  voids  and  cleaved  surfaces.  Scanning  electron  microscopy  and  energy  dispersive 
spectroscopy  reveal  uniform  distribution  of  PTFE,  providing  nanoscale  mixing  within  particles.  The 
combustion  enthalpy  was  found  to  be  20.2  kJg-1,  though  a  slight  decrease  (0.8  kJg-1)  results  from 
extended  high  energy  milling  due  to  a-AlF3  formation.  For  high  energy  mechanically  activated  particles, 
differential  scanning  calorimetry  in  argon  shows  a  strong,  exothermic  pre-ignition  reaction  that  onsets 
near  440  °C  and  a  second,  more  dominant  exotherm  that  onsets  around  510  °C.  Scans  in  02-Ar  indicate 
that,  unlike  physical  mixtures,  more  complete  reaction  occurs  at  higher  heating  rates  and  the  reaction 
onset  is  drastically  reduced  (approx.  440  °C).  Simple  flame  tests  reveal  that  these  altered  Al- 
polytetrafluoroethylene  particles  light  readily  unlike  micrometer- sized  aluminum.  Safety  testing  also 
shows  these  particles  have  high  electrostatic  discharge  (89.9-108  mJ),  impact  (>213  cm),  and  friction 
(>360  N)  ignition  thresholds.  These  particles  may  be  useful  for  reactive  liners,  thermobaric  explosives, 
and  pyrolants.  In  particular,  the  altered  reactivity,  large  particle  size  and  relatively  low  specific  surface 
area  of  these  fuel  rich  particles  make  them  an  interesting  replacement  for  aluminum  in  solid  propellants. 
See  DOI:  10. 1002/prep. 20 1200 102  for  more  details.  Fluorinated  graphite,  poly(carbon  monofluoride) 
(PMF),  was  also  considered  as  an  inclusion  material  because  of  its  vastly  different  mechanical  properties. 
See  DOI:  10. 1002/prep. 20 1200202  for  more  details. 

V.  Investigation  of  Solid  Oxidizer  and  Gaseous  Fuel  Combustion  Performance  Using 
an  Elevated  Pressure  Counterflow  Experiment  for  Reverse  Hybrid  Rocket  Engine  (Penn 
State  including  Grant  Risha,  Penn  State- Altoona,  and  Greg  Young,  Naval  Surface  Warfare 
Center,  Indian  Head) 

Pressurized  counterflow  burner  and  static-fired  motor  studies  were  conducted  to  explore  the 
possibility  of  a  reverse  hybrid  system,  having  a  solid  oxidizer  and  gaseous  fuel.  Theoretical  performance 
analysis  indicates  such  a  system  may  yield  specific  impulse  and  density  specific  impulse  similar  to 
composite  solid  propellants,  yet  with  the  added  capability  to  throttle,  shut  down,  and  restart  firings. 
Pressurized  counterflow  flame  studies,  conducted  using  pressed  ammonium  perchlorate  (AP)  pellets  and 
gaseous  ethylene,  show  three  pressure  dependent  combustion  regimes.  These  results,  along  with 
published  results  from  Hightower  and  Price  [188],  Boggs  [189],  and  Atwood  et  al.  [190]  are  presented  in 
Figure  65. 
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Figure  65  Counterflow  burning  rates  as  a  function  of  pressure,  compared  to  published  data 
(Hightower  and  Price  [188],  Boggs  [189],  Atwood  et  al.  [190]). 

Ammonium  perchlorate  decomposition,  for  pressures  below  1  MPa,  is  controlled  by  heat  transfer  from  the 
resulting  fuel/oxidizer  diffusion  flame,  which  forms  between  the  fuel  and  decomposition  products  of  AP. 
In  this  low  pressure  regime,  the  AP  burning  rate  is  found  to  increase  with  flame  strain  rate  and  pressure, 
yielding  measured  values  between  0.1  to  0.5  mm/s.  As  pressure  increases,  the  monopropellant  flame 
moves  closer  to  the  oxidizer  surface  until  the  pressure  reaches  the  self-decomposition  limit  near  3  MPa,  at 
which  point  the  monopropellant  flame  becomes  nearly  independent  of  the  diffusion  flame.  Further 
increasing  the  pressure  yields  burning  rates  between  0.4  to  0.7  cm/s,  which  are  consistent  with  the 
literature.  Variation  of  flame  strain  rate  under  these  conditions  has  little  or  no  influence  on  the  AP 
burning  rate  for  the  range  of  flow  conditions  tested.  Similar  studies  conducted  with  methane  suggest 
burning  rates  are  not  affected  by  fuel  type.  Lab-scale  static  motor  firings  were  conducted  to  examine 
ignition,  variation  of  fuel  flow  rate  and  initial  motor  pressure,  and  system  performance.  Results  indicate 
that  successful  motor  operation  requires  initial  pressures  capable  of  boosting  the  system  into  the  higher 
burning  rate  regimes  of  AP.  Ignition  at  low  pressures  results  in  long  startup  periods  where  the  motor 
operates  as  a  gas  generator. 

VI.  Hierarchical  Theoretical  Methods  for  Understanding  and  Predicting  Anisotropic 
Thermal  Transport  and  Energy  Release  in  Rocket  Propellant  Formulations  (Thomas  D. 
Sewell  (PI),  Chemistry,  U.  of  Missouri-Columbia) 

There  is  a  strong  push  to  explore  possibilities  for  new  explosives  and  propellants  based  on  novel 
combinations  of  traditional  polymer  fuels,  ionic  salt  oxidizers,  metals,  metal  oxides,  and  CHNO  energetic 
materials  and  less  traditional  materials  such  as  nanoparticles,  thermites,  graphene,  carbon  nanotubes 
(CNTs),  and  boron,  among  others.  While  equilibrium  thermodynamics  calculations  yield  estimates  of 
overall  theoretical  energy  release  for  a  given  propellant  formulation  with  arbitrary  combinations  of  both 
traditional  and  exotic  ingredients,  they  do  not  account  for  non-equilibrium  spatio-temporal  processes  of 
transport  and  chemistry  on  “sub-continuum”  spatial  scales  and  mesoscopic  interactions  among  the 
propellant  ingredients  that  determine  the  energy  release  rates  and  performance.  When  gross  theoretical 
considerations  based  on  simple  thermodynamic  equilibrium  are  used  to  predict  specific  impulse  and 
stability  criteria  for  a  given  propellant  formulation,  those  estimates  often  are  poor  and  hence  deemed  not 
predictive  when  compared  with  results  from  experimental  measurement  and  testing. 

One  approach  for  achieving  theoretically  based  performance  limits  is  to  more  carefully  prescribe 
the  interactions  of  the  propellant  ingredients  in  intimate  contact  and  the  nature  of  mixing  and  chemistry  of 
reactants  and  products  during  microstructural  decomposition  and  burning.  Improved  and  systematic 
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descriptions  that  consider  detailed  non-equilibrium  processes  (when  viewed  at  the  macroscale  or  system 
level  scale)  can  lead  to  better  predictions  and  enable  the  generation  of  new  propellant  designs.  New 
designs  aimed  at  enhanced  performance  will  use  smaller  and  more  precisely  engineered  particles, 
approaching  the  nanometer  scale,  to  increase  the  proximity  between  interacting  propellant  ingredients; 
use  orientation,  placement,  and  the  intrinsic  anisotropy  of  material  behaviors  of  the  individual  ingredients 
as  a  path  toward  thermal  management  of  the  system  in  the  pre-  or  early-reaction  stages;  and  use 
engineered  microstructures  on  the  representative  volume  element  (RVE)  scale  and  above  to  achieve 
desired  mixing  and  reactive  flow  properties  on  the  mm  scale. 

An  attractive  approach  for  developing  new  designs  is  to  use  self-  and  direct-assembly  to  generate 
highly  controlled  nano-  and  microstructure  that  will  typically  yield  close  intimacy  but  with  additional 
degrees  of  freedom,  corresponding  to  placement  and  orientation  of  ingredients,  that  will  enable  precise 
tailoring  of  properties.  Such  materials  are  usually  layered  or  tiled,  with  well-defined  orientable  crystalline 
or  quasi-crystalline  properties.  In  general,  such  directed  assemblies  will  exhibit  anisotropy  that  extends 
beyond  that  associated  with  the  individual  ingredients.  This  allows  one  to  exploit  the  intrinsic  anisotropy 
in  the  mechanical  and  electromagnetic  properties  of  individual  propellant  ingredients  and  moreover  to 
impose  hierarchical  design  at  the  nano-  and  microscales,  a  feature  we  refer  to  as  hierarchical  anisotropy. 
Together,  these  should  make  it  possible  to  cause  directionality  in  the  overall  reaction  and  energy  release 
profile  in  propellants.  Some  possible  methods  of  achieving  this  that  we  will  study  include  use  of 
anisotropic  assemblies  such  as  bundles  of  reacting  linear  chains  of  nanoparticles,  anisotropic  wave 
propagation  engineered  in  via  anisotropic  crystalline  structures,  or  anisotropic  combustion  engineered 
using  highly  anisotropic  thermal  conductors  such  as  graphene  sheets. 

The  efficacy  of  new  propellants  based  on  hierarchically  anisotropic  materials  depends  on 
structures  in  the  unreacted  material,  but  also  on  the  influence  of  these  structures  in  the  decomposing  and 
burning  propellant.  Therefore,  it  is  useful  to  think  of  anisotropy  at  three  levels: 

1)  Intrinsic  anisotropy  at  the  molecular  up  to  the  continuum  microscale  for  pure  constituents 

2)  Manufactured  nano-  and  microscale  anisotropy  that  is  induced  by  manufacture  specifications  of 
the  composition  that  identifies  the  material  selection,  gross  stoichiometry,  morphology,  and 
constituent  placement  and  orientation 

3)  Mesoscale  anisotropy  persistence  during  the  physico-chemical  structural  decomposition,  mixing, 
and  reactive  processes,  templated  by  item  2),  that  results  in  locally  anisotropic  effects  that  appear 
homogenized  at  the  engineering  device-level  scale  that  ultimately  defines  the  performance 
characteristics  of  the  propellant. 

The  goal  of  our  effort  is  a  focused,  integrated  program  of  fundamental  theoretical  and 
computational  research  designed  to  yield  understanding  and  predictive  capability  regarding  anisotropic 
thermal  transport  and  energy  release  in  advanced  rocket  propellants.  The  proposed  research  will 
contribute  to  a  practical  capability  for  the  a  priori  design  of  smart,  functional  propellant  materials  based 
on  intelligently  designed  polymer  nanocomposite  formulations  augmented  by  non-traditional  additives  or 
passivation  agents,  possibly  coupled  to  external  fields. 

Our  approach  is  to  develop  a  hierarchical  theoretical  framework  that  combines,  with  as  much 
rigor  as  possible: 

•  Fundamental  information  from  atomic-scale  simulations 

•  Detailed  continuum-based  mesoscopic  models  of  interfaces  and  interphases,  phase 
transformations,  multi-phase  transport,  and  chemistry 

•  Explicit  microstructure-resolved  RVE  scale  and  larger  simulations  that  generate  homogenized 
results  induced  by  mesoscale  constituents 

•  Development  of  approaches  that  define  engineering-scale  models  for  propellant  combustion 
uniquely  defined  by  the  hierarchical  structures  from  the  molecular,  micro-,  and  mesoscopic 
models.  These  engineering-scale  models  are  the  ones  that  can  be  used  to  make  system 
performance  estimates. 
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Atomic-scale  methods  -  molecular  dynamics  (MD)  and  electronic  structure  theory  -  are  being 
used  to  predict  the  details  of: 

•  Thermo-mechanical  and  transport  properties 

•  Thermodynamic  phase  diagrams 

•  Interfacial  structure  and  processes 

•  Ignition  chemistry  in  multiple  phases  and  mixing  regimes. 

The  simulations  are  being  done  for  conditions  of  temperature,  stress-/strain-state,  and  tensoral  strain  rate 
relevant  to  propellant  ignition  and  initiation  and  are  being  provided  to  the  mesoscale  and  continuum  scale 
theorists  for  direct  use  in  their  models/predictions. 

Mesoscopic  continuum-based  models  are  being  developed  for  chemically  reactive  materials.  The 
general  approach  taken  follows  standard  methodologies  used  in  the  combustion  community,  but  includes 
the  complexities  associated  with  anisotropy  and  heterogeneity  that  must  be  considered  to  describe  nano- 
and  microscale  physical  processes  and  their  effects  on  chemistry.  Those  include  consideration  of  phase 
changes  including  solid-state  polymorphism  and  transitions  among  solid,  liquid,  and  gas;  and  tracking  of 
key  species  diagnostic  of  or  critical  to  the  overall  reaction  kinetics.  Multiphase  diffusion  of  species  and 
diffusion  of  oxidizing  and  metal  ions  in  oxide  melt  solutions  (often  thought  to  be  involved  in  metal, 
intermetallic,  and  thermitic  reactions)  are  being  included.  The  mesoscopic  models  will  predict: 

•  Multiphase  decomposition  and  reaction  during  chemical  energy  release  at  initially  separate 
diffusive/reactive  interfaces 

•  Volumetric  reaction  and  gaseous  combustion 

•  Formalism  and  constitutive  theory  used  to  describe  reactants  and  their  products  that  reflect  the 
underlying  material  mechanics  and  physical  chemistry  of  the  constituents 

•  Emission  outputs  (spectra)  for  “tagged”  species  in  the  continuum  simulations  that  can  be  used  to 
generate  synthetic  spectra  for  comparison  to  experiment. 

Thus,  we  will  generate  highly  constrained  model  simulations  of  multiphase  decomposition  and 
combustion  events  for  which  model  spectra  output  can  be  compared  directly  to  both  experiments  and  MD 
simulations  to  develop  temporally  and  spatially  accurate  reaction  kinetics. 

On  the  basis  of  the  mesoscopic  models  we  are  developing  RVE-scale  methods  for  simulating 
experimentally  relevant  configurations  such  as: 

•  1-,  2-  and  3-D  laminates  and  plies 

•  Spherical  core/shell  particles 

•  Both  periodic  and  stochastic  microstructure  with  designated  distributions  of  materials. 

Ultimately  we  will  simulate  a  piece  of  material  sufficiently  large  to  yield,  in  a  proper  homogenization 
over  the  ’’microscopic”  ensemble,  statistical  distributions  of  properties  that  can  be  used  in  macro-scale 
models. 

During  the  third  year  of  the  overall  project,  considerable  effort  was  devoted  to  developing  the 
fundamental  methods  to  be  applied  in  the  successive  years  of  the  project.  The  work  is  essentially  on  track 
compared  to  what  was  proposed  originally.  A  brief  summary  from  each  participating  institution  is 
provided  below. 


A.  Development  and  implementation  of  atomic-scale  methods  for  computing  thermal- 
mechanical-chemical  properties  (Thomas  D.  Sewell  and  Donald  L.  Thompson,  U.  of  Missouri- 
Columbia) 

Efforts  in  the  third  year  were  focused  on  the  development  and  implementation  of  atomic-scale 
methods  for  computing  thermal-mechanical-chemical  properties  that  will  be  provided  directly  to  the 
mesoscopic  (Matalon  and  Stewart  at  UIUC)  and  macroscopic  (Ortiz  at  Caltech)  modelers  on  the  project. 
The  materials,  properties,  and  thermodynamic  states  were  chosen  in  cooperation  with  those  modelers.  The 
major  areas  of  emphasis  at  the  University  of  Missouri-Columbia  during  the  this  year  were: 

•  Anisotropic  thermal  conductivity  in  triclinic  crystalline  compounds 
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•  Second-order  elastic  tensor  for  triclinic  crystalline  compounds  for  finite  temperatures  and 
pressures 

•  Calculations  of  plasticity  in  the  basal  plane  of  the  crystalline  energetic  material  TATB 

•  Melting  point  of  anisotropic  crystals,  both  at  normal  and  elevated  pressure 

•  Liquid  state  transport  coefficients  as  functions  of  temperature  and  density 

•  Force  field  development  for  TATB  and  FOX-7 

•  Electronic  structure  predictions  of  the  transformation  from  acid-base  pairs  (e.g.,  nitric  acid  and 
ammonia)  to  ion  pairs  (e.g.,  NH4+  and  N03“),  that  is,  proton  transfer,  in  small  clusters  of 
ammonium  nitrate  (AN)n  and  ammonium  perchlorate  (AP)n  containing  n  >  1  formula  units  of  AN 
or  AP. 

•  Electronic  structure  calculations  of  the  oxidation  of  A1  by  C02 

Molecular  dynamics  (MD)  and  molecular  mechanics  (MM)  simulations  were  used  to  predict  the 
thermal  and  mechanical  properties  of  the  insensitive  energetic  material  1, 3, 5-triamino-2, 4,6- 
trinitrobenzene  (TATB)  in  both  the  crystalline  and  liquid  states.  Although  it  is  not  used  as  a  propellant 
ingredient,  TATB  crystallizes  in  a  triclinic  space  group  (the  lowest  symmetry  class)  and,  therefore,  codes 
and  analysis  tools  developed  for  TATB  should  be  applicable  to  all  other  molecular  energetic  materials. 
Progress  was  also  made  toward  the  development  of  a  flexible-molecule  force  field  for  1,1 -diamino-2, 2- 
dinitroethene  (FOX-7),  which  is  of  potential  interest  as  a  propellant  additive. 

Proton  transfer  is  often  thought  to  be  an  early  reaction  in  the  decomposition  of  inorganic  salts 
such  as  ammonium  nitrate  (AN)  and  ammonium  perchlorate  (AP),  which  are  commonly  used  oxidizers. 
As  an  initial  stage  of  developing  chemical  decomposition  schemes  for  AN  and  AP,  electronic  structure 
(quantum  chemistry)  calculations  were  performed  to  study  the  transformation  from  the  stable  acid-base 
pair  for  isolated  formula  units  to  stable  ion  pairs,  as  described  in  the  bulleted  list  just  above.  In  addition  to 
calculations  for  pure  (AN)n  and  (AP)n,  the  effects  of  NH3,  A1H3,  and  BH3  molecules  on  cluster  acid- 
base/ion-pair  stability  were  studied. 

B.  Development  of  Mesoscale  Simulations  and  Theory  for  Multi-Material  and 
Polycrystalline  and  Microstructured  Materials  (Michael  Ortiz,  Caltech) 

The  work  tasked  to  Caltech  is  concerned  with  the  development  of  mesoscale  simulations  and 
theory  for  multi-material  and  polycrystalline  and  microstructured  materials.  Specific  subtasks  are  the 
generation  of  useable,  verified  and  validated  material  descriptions,  mixtures  and  interfacial  material 
interactions  that  can  be  directly  incorporated  into  any  larger  continuum-based  computational  framework. 
A  primary  focus  of  the  work  specifically  concerns  the  development  of  hierarchical  theoretical  methods  for 
understanding  and  predicting  anisotropic  thermal  transport  and  energy  release  in  rocket  propellant 
formulations.  An  even  more  specific  focus  concerns  the  development  of  understanding  and  predictive 
capability  regarding  anisotropic  thermal  transport  and  energy  release  in  advanced  rocket  propellants, 
leading  to  a  practical  capability  for  a  priori  propellant  design.  Propellants  have  traditionally  been 
modeled  using  hydrocodes  that  neglect  entirely  the  strength  of  the  material.  However,  recent  work 
(Kroonblawd  &  Sewell,  J.  Chem.  Phys .,  2013)  has  shown  that  the  thermal  conductivity  of  energetic 
materials  such  as  TATB  can  be  extremely  anisotropic,  with  conductivity  contrasts  of  orders  of  magnitude 
depending  on  transport  direction.  It  seems  likely  that  this  extreme  anisotropy  should  in  turn  result  in  a 
strong  directional  dependence  of  reaction-front  speeds,  a  phenomenon  that  appears  to  be  poorly 
understood  at  present.  This  acute  sensitivity  to  local  anisotropy  in  turn  raises  opportunities  for  the  optimal 
design  of  microstructures  resulting  in  enhanced  engineering  properties  of  propellants,  such  as  specific 
impulse,  a  design  tradeoff  that  also  appears  to  be  unexplored  at  present. 

Work  at  Caltech  to  date  has  progressed  on  two  different  fronts,  described  next  in  turn. 
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Figure  66  The  paradigm  shift:  From  Newtonian  dynamics  to  diffusional  transport  (heat  and  mass) 

Time  step  limited  by  diffusional  time  scale 

Fundamental  understanding  of  the  effect  of  anisotropy  on  reaction-front  speeds.  The  essential  difficulty 
in  developing  such  understanding  is  the  multiscale  nature  of  the  phenomenon,  which  includes  both  atomic 
level  rate-limiting  processes,  such  as  thermal  vibrations  and  bond  breaking,  collisions,  coupled  to  slow 
processes  such  as  heat  and  mass  transfer,  viscosity  and  complex  reactions  paths.  Indeed,  the  challenge  of 
performing  atomistic  full-chemistry  simulations  over  time-scales  relevant  to  macroscopic  chemical 
processes  and  diffusion-mediated  phenomena  constitutes  a  long-standing— yet  stubbornly  elusive — goal 
in  computational  science.  For  propellants,  the  time-scale  gap  is  staggering,  from  thermal  vibrations 
(femptosecond)  to  device  operation  (seconds).  The  spatial-scale  gap  is  equally  staggering,  from  atomic 
lattice  scale  (Angstroms)  to  device  dimensions  (m).  And  yet,  neither  of  the  limiting  scales  is  sufficient  to 
provide  a  good-enough  handle  for  material  and  device  design:  On  one  had  we  wish  atomistic  realism 
across  within  reaction  zone  in  order  to  account  for  complex  difficult-to-predict  reaction  paths;  on  the 
other  hand  we  wish  to  make  contact  with  engineering  systems  and  data.  The  fundamental  question  is, 
therefore:  How  to  affect  space-time  coarse-graining  (atoms  to  device)  without  introducing  spurious 
physics  and  without  essential  loss  of  information?  In  effect,  we  wish  for  atomistic  realism  without  the 
need  to  resolve  the  length  scale  of  the  crystal  lattice  and  the  time  scale  of  thermal  vibrations.  We  have 
addressed  this  challenge  by  jettisoning  ab  initio  molecular  dynamics  if  favor  of  non-equilibrium  statistical 
mechanics.  In  this  paradigm  shift,  phase  space  is  replaced  by  atomic-level  temperatures  and  molar- 
fractions,  Newtonian  dynamics  by  atomic-level  heat  and  mass  diffusive  transport,  the  required  spatial 
resolution  is  controlled  by  temperature  and  concentration  gradients  and  the  required  temporal  resolution  is 
controlled  by  diffusional  transients.  The  entire  theory  is  derived  consistently  from:  Jayne’s  maximum 
entropy  principle,  which  enables  the  definition  of  atomic-level  free  entropies  with  non-uniform 
temperatures  and  molar  fractions;  and  from  a  reformulation  of  Onsager  kinetics  at  the  atomistic  level.  The 
approach  supplies  atomistic  realism  and,  in  particular,  it  takes  as  input  an  arbitrary  (empirical)  interatomic 
potential.  In  addition,  empirical  kinetic  equations  may  also  be  supplied  at  the  atomistic  level.  The  spatial 
resolution  of  the  model  is  additionally  adapted  according  to  local  resolution  needs,  from  atomistic  to 
continuum,  using  the  quasicontinuum  method  (E.B.  Tadmor  et  al. ,  Phil  Mag.  A ,  73(6):  1529-1563, 
1996).  We  are  in  the  process  of  implementing  the  approach  and  applying  it  to  a  simple  validation/test 
problem:  the  combustion  of  graphite.  Graphite  supplies  a  good  test  material  as  it  exhibits  two-order-of- 
magnitude  contrast  in  the  components  of  its  thermal  conductivity  tensor.  We  specifically  simulate  a  slab 
of  graphite  in  contact  with  an  oxygen  reservoir,  and  track  the  combustion  reaction  as  it  progresses  through 
the  slab,  the  evolution  of  the  temperature  field,  the  elasticity  of  the  unreacted  graphite  and  the 
hydrodynamics  of  the  reaction  products.  In  order  to  account  for  chemical  reactions  we  use  a  ReaxFF  C-O 
interatomic  potential  available  in  LAMMS.  The  aim  is  to  ascertain  the  fidelity  of  predictions  of  the 
dependence  of  the  reaction-front  speed  on  the  orientation  of  the  graphite  lattice. 
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Figure  67  Propellant  Microstruture 


Development  of  methods  to  optimize  propellant  microstructures.  The  aim  of  this  effort  build  on 
fundamental  information  from  atomic-scale  calculations,  including  orientation-dependent  reaction-front 
speeds,  in  order  to  optimize  parameterized  propellant  microstructures  of  interest,  including  composites 
comprising  a  Polymer  binder,  oxidizer  granules  and  additives  (e.g.,  NH4CIO4  Composite  Propellant, 

APCP).  Key  design  parameters  to  be  optimized  include  volume  fractions,  textures  and  grain/inclusion 
morphologies.  We  connect  micromechanical  properties  to  overall  propellant  performance  through  front 
tracking  simulations  of  the  passage  of  reaction  fronts  through  representative  volume  elements,  thus 
determining  effective/macroscopic  burning  rates.  Such  effective  burning  rates  are  subsequently  optimized 
for  maximum  specific  impulse.  To  date,  we  have  completed  the  development  and  testing  of  a  ray-tracing 
solver  with  capability  for  tracking  the  passage  of  reaction  fronts  through  arbitrary  locally-anisotropic 
spatial  microstructures.  At  present,  the  local  reaction-front  speed  depends  on  the  local  anisotropic  thermal 
conductivities  through  a  simple  one-step  reaction  model,  but  the  computational  capability  is  sufficiently 
general  to  integrate  micromechanical  models  as  they  become  available.  The  front  tracking  calculations  are 
fast  enough  as  to  make  microstructural  optimization  feasible. 

C.  Continuum  Modeling  and  Simulation  of  Microscopic,  Multi-phase,  Reactive 
Processes  of  Reactants  at  the  Nano-  and  Microscale  (D.  Scott  Stewart  andMoshe  Matalon, 

UIUC) 

The  efforts  of  Professors  D.  Scott  Stewart  and  Moshe  Matalon,  and  their  students  have  been  directed  at 
continuum-based  modeling  and  simulation  of  multifunctional  energetic  materials  of  interest  to  the  US  Air  Force 
Office  of  Scientific  Research,  for  advanced  propulsion.  For  the  contract  year,  the  main  effort  has  been  modeling 
of  condensed  phase  reaction/diffusion  of  mobile  species  with  applications  to  metal  combustion.  And,  specifically, 
the  Stewart/Matalon  team  developed: 

1)  A  model  of  generic  solid/liquid  condensed  phase  combustion  of  separated  reactants  (reaction/diffusion 
flames)  dubbed  the  ’’Gibbs  formulation”.  This  included  modeling  of  both  energetic  solids  and  liquids.  The 
work  on  the  multicomponent  formulation  is  being  prepared  for  publication  [191]. 

2)  Developed  an  analysis  of  a  fundamental  small-scale  three  components  model  (that  is  applicable  to  metal 
combustion  with  oxidizers).  This  model  has  been  articulated  for  Aluminum  (Al),  oxidizers  and  binders  in 
particular  in  a  generic  counterflow  flame  configuration.  This  work  is  currently  being  prepared  for 
publication  [192]. 

3)  Started  the  development  of  a  theory  of  a  condensed  phase  homogenization/averaging. 

Modeling  of  condensed  phase  reaction/diffusion  of  mobile  species,  with  applications  to  metal  combustion 
(Matalon/Stewart/Bdzil).  Solid  propellants  are  mixtures  of  metal  particles,  binders  and  various  oxidizers 
that  are  combined  in  a  mixture  that  is  molded  to  manufacture  large  propellant  grain  configurations. 
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Substantial  decomposition  of  the  reactive  materials  takes  place  as  the  gross  gas/condensed  phase  interface 
(surface)  progresses  through  the  large  propellant  grain.  Exothermic  reaction  takes  place  in  the  condensed 
phase  in  the  sub-surface  regions.  Moreover,  the  extreme  temperature  gradients  and  the  associated  thermal 
expansion  cause  the  condensed  phase  materials  to  swell  and  set  up  flows  around  reactants  that  are  initially 
separated.  At  the  interfaces  of  reactants  (metals  and  binders)  and  oxidizers  (like  Ammonium  Perchlorate 
or  other  novel  oxidizers)  diffusion  flames  in  the  condensed  materials  can  be  established. 

Initial  considerations  included  looking  at  generic  configurations  that  are  representative  of  actual 
or  realizable  micro  structures.  Figure  68  shows  a  representative  microstructure  of  a  composite  solid 
propellant.  Figure  69  shows  a  characteristic  schematic  (based  on  a  real  microstructure)  that  indicates  the 
generic  and  ubiquitous  nature  of  a  counterflow  configuration  that  has  two  opposed  flows  of  reactants  (say 
Fuel  and  Oxidizer)  in  the  condensed  phase  as  one  nears  the  condensed/gas  surface.  Figure  70  shows  an 
inset  that  highlights  the  counterflow  configuration  our  team  analyzed  in  this  past  year 
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Figure  68  A  representative  microstructure  of  a  composite  solid  propellant. 


Decomposition  and  reaction  in  the  condensed  phases 
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There  are  relative  motions  in  the  molten  layer  due  to  thermal  expansion  resulting 
from  the  large  temperature  change  (up  to  3000  K),  with  velocities  comparable 
to  the  receding  front,  which  is  moving  at  ^  1  cm/s. 


Figure  69  Schematic  of  a  temporary  state  during  the  burning  of  a  solid  propellant  (based  on  the 
microstructure  shown  in  Fig.  68).  The  dashed  corrugated  curve  on  top  illustrates  the  interface 

advancing  through  the  propellant  grain. 
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Multi-component,  model  for  separated  reactants 
Condensed  phase  diffusion  flames 


Figure  70  The  inset  on  the  right  shows  the  relevance  of  the  counterflow  configuration  studied  by  the 
Matalon/Stewart  team  in  the  contract  year  to  the  large-scale  propellant  burning. 

Our  team  analyzed  the  counter  flow  flame  configuration  that  has  the  following  model  elements: 

•  General  equations  for  a  multi-component  (condensed)  reactive  material 

•  A  Stefan-Maxwell  diffusion  model  that  leads  to  a  generalized  Fick's  diffusion  model 

•  Temperature  and  concentration  dependent  diffusion 

•  A  generic  three  component  counter  flow  flame  configuration  with  2  reactants  that  combine  to 
generate  products 

•  Arrhenius  reaction  rate  kinetics 

•  Asymptotic  analysis  for  large  Damkohler  number  and  activation  energy,  and  explicit  flame 
structure  formulas 

•  Accompanying  numerical  solutions  of  the  model  equations  for  verification  and  generalization 
purposes 

The  work  represented  by  the  above  list  is  being  prepared  for  publication  [192,193]. 

Figure  71  is  a  sketch  of  the  flame  structure  that  illustrates  the  unique  characteristics  of  diffusion 
flames  in  condensed  reactants.  The  striking  element  is  related  to  the  very  small  diffusion  mobility  of  the 
condensed  species  which,  in  contract  to  gaseous  flames,  can  only  diffuse  towards  the  reaction  zone  is  a 
very  thin  mixing-zone  where  the  temperature  is  sufficiently  high.  The  mixing  layer  is  seen  to  be 
significantly  narrower  than  the  heat  conduction  layers  that  are  due  to  the  0(1)  thermal  diffusivities,  and 
extend  to  much  larger  distances.  In  addition  this  type  of  flame  leads  to  the  possibility  of  extinction  when 
the  strain  rate  exceeds  a  critical  value  (or  the  properly  defined  Damkohler  number  D  is  reduced 
significantly,  see  Figure  72).  This  picture  indicates  the  basic  result  that  is  being  prepared  for  publications 
[192]  &  [193]. 
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Flame  structure 


reaction  zone 


0(0- l) 


heat  conduction  region 


Figure  71  Sketch  of  the  generic  condensed  phase  diffusion  flame  structure  based  on  recent 

asymptotic  results  by  the  UIUC  team. 


State  diagram  -  burning/non-burning  states 
ignition/extinction  phenomena 


Damkohler  number  JD  ^  residence  time 
controlled  by  the  strain  rate  e. 

Figure  72  Typical  response  curve  of  flame  temperature  vs  Damkohler  number,  illustrating  the 
possibility  of  extinction  states  at  sufficiently  high  strain  rates  (small  D),  based  on  our  recent 

asymptotic  and  numerical  results,  [192]. 

Development  of  a  continuum  model  for  an  anisotropic  energetic  material  that  undergoes  phase  change 
and  chemical  reaction  subject  to  thermal  loading  and  shear  stress.  As  part  of  the  modeling  effort,  the 
general  multi-component  formulation  (dubbed  the  ’’Gibbs  formulation”)  has  been  extended  to  develop  a 
model  for  energetic  solids  that  undergo  decomposition.  Some  candidates  are  TATB  crystals  and  Fox  7 
(DADNE).  This  work  is  being  done  in  collaboration  with  the  U  Missouri  team.  Figure  73  shows  an 
overview  sketch  that  represents  this  effort. 
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B)  Development  of  a  continuum  model  for  an  anisotropic 
energetic  material  that  undergoes  phase  change  and 
chemical  reaction  subject  to  thermal  loading  and  shear 
stress. 

Mixture  of  a  small-strain  (anisotropic)  elastic  solid  and 
compressible  melt  liquids  and  products 

A  — >•  B  — >•  C 

solid  liquid  Products  (I) 

Yi  Y2  Y, 

Figure  73  The  UIUC  team  is  developing  a  continuum-based  model  informed  by  atomistics  for 

anisotropic  energetic  solids. 
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